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METHOD FOR OPTIMIZING A PLANT WITH 
MULTIPLE INPUTS 

CROSS-REFERENCE TO RELATED 
APPLICATIONS 

[0001] This application is a Continuation of U.S. patent 
application Ser. No. 09/224,648, filed Dec. 31, 1998, and 
entitled "AMETHOD FOR OPTIMIZING A PLANT WITH 
MULTIPLE INPUTS," (Atty. Dkt. No, PAVI-24,541), which 
is a Continuation-in-Part Application of U.S. patent appli- 
cation Ser. No. 09/167,504, filed Oct. 6, 1998, and entitled 
"A METHOD FOR ONE-LINE OPTIMIZATION OF A 
PLANT," (Atty Dkt No. PAVI-24100), which is a Continu- 
ation-in-Part Application of U.S. patent application Ser. No. 
08/943,489, filed Oct. 3, 1997, and entitled "A METHOD 
FOR STEADY-STATE IDENTIFICATION BASED UPON 
IDENTIFIED DYNAMICS," which is a Continuation-in- 
Part of U.S. patent Ser. No. 08/643,464, filed May 6, 1996, 
and entitled "Method and Apparatus for Modeling Dynamic 
and Steady-State Processors for Prediction, Control, and 
Optimization." 

TECHNICAL FIELD OF THE INVENTION 

[0002] The present invention pertains in general to neural 
network based control systems and, more particularly, to 
on-line optimization thereof 

BACKGROUND OF THE INVENTION 

[0003] Process models that are utilized for prediction, 
control and optimization can be divided into two general 
categories, steady-state models and dynamic models. In each 
case the model is a mathematical construct that characterizes 
the process, and process measurements are utilized to 
parameterize or fit the model so that it replicates the behav- 
ior of the process. The mathematical model can then be 
implemented in a simulator for prediction or inverted by an 
optimization algorithm for control or optimization. 

[0004] Steady-state or static models are utilized in modem 
process control systems that usually store a great deal of 
data, this data typically containing steady-state information 
at many different operating conditions. The steady-state 
information is utilized to train a non-linear model wherein 
the process input variables are represented by the vector U 
that is processed through the model to output the dependent 
variable Y The non-linear model is a steady-state phenom- 
enological or empirical model developed utilizing several 
ordered pairs (U h YJ of data from different measured steady 
states. If a model is represented as: 



Y«?(u,Y) 



(001) 



[0005] where P is some parameterization, then the 
steady-state modeling procedure can be presented as: 



(U,7)-P 



(002) 



[0006] where U and Y are vectors containing the U ( , 
Y| ordered pair elements. Given the model P, then the 
steady-slate process gain can be calculated as: 



A/W Y) 



(003) 



[0007] The steady-state model therefore represents the 
process measurements that are taken when the system is in 
a "static** mode. These measurements do not account for the 
perturbations that exist when changing from one steady-state 
condition to another steady-state condition. This is referred 
to as the dynamic part of a model. 

[0008] A dynamic model is typically a linear model and is 
obtained from process measurements which are not steady- 
state measurements; rather, these are the data obtained when 
the process is moved from one steady-state condition to 
another steady-state condition. This procedure is where a 
process input or manipulated variable u(t) is input to a 
process with a process output or controlled variable y(t) 
being output and measured. Again, ordered pairs of mea- 
sured data (u(I), y(I)) can be utilized to parameterize a 
phenomenological or empirical model, this time the data 
coming from non-steady-state operation. The dynamic 
model is represented as: 

y(t)=p(u(t)#(t)) (004) 

[0009] where p is some parameterization. Then the 
dynamic modeling procedure can be represented as: 

(u*,7)-»P (005) 

[0010] Where u and y are vectors containing the (u(I),y(I)) 
ordered pair elements. Given the model p, then the steady- 
state gain of a dynamic model can be calculated as: 



k = 



(006) 



[0011] Unfortunately, almost always the dynamic gain k 
does not equal the steady-state gain K, since the steady-state 
gain is modeled on a much larger set of data, whereas the 
dynamic gain is defined around a set of operating conditions 
wherein an existing set of operating conditions are mildly 
perturbed. This results in a shortage of sufficient non-linear 
information in the dynamic data set in which non-linear 
information is contained within the static model. Therefore, 
the gain of the system may not be adequately modeled for an 
existing set of steady-state operating conditions. Thus, when 
considering two independent models, one for the steady- 
state model and one for the dynamic model, there is a 
mis-match between the gains of the two models when used 
for prediction, control and optimization. The reason for this 
mis-match are that the steady-state model is non-linear and 
the dynamic model is linear, such that the gain of the 
steady-state model changes depending on the process oper- 
ating point, with the gain of the linear model being fixed. 
Also, the data utilized to parameterize the dynamic model do 
not represent the complete operating range of the process, 
i.e., the dynamic data is only valid in a narrow region. 
Further, the dynamic model represents the acceleration 
properties of the process (like inertia) whereas the steady- 
state model represents the tradeoffs that determine the 
process final resting value (similar to the tradeoff between 
gravity and drag that determines terminal velocity in free 
fall). 

[0012] One technique for combining non-linear static 
models and linear dynamic models is referred to as the 
Hammerstein model. The Hammerstein model is basically 
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an input-output representation that is decomposed into two 
coupled parts. This utilizes a set of intermediate variables 
that are determined by the static models which are then 
utilized to construct the dynamic model. These two models 
are not independent and are relatively complex to create. 

[0013] Plants have been modeled utilizing the various 
non-linear networks. One type of network that has been 
utilized in the past is a neural network. These neural 
networks typically comprise a plurality of inputs which are 
mapped through a stored representation of the plant to yield 
on the output thereof predicted outputs. These predicted 
outputs can be any output of the plant. The stored represen- 
tation within the plant is typically determined through a 
training operation. 

[0014] During the training of a neural network, the neural 
network is presented with a set of training data. This training 
data typically comprises historical data taken from a plant. 
This historical data is comprised of actual input data and 
actual output data, which output data is referred to as the 
target data. During training, the actual input data is presented 
to the network with the target data also presented to the 
network, and then the network trained to reduce the error 
between the predicted output from the network and the 
actual target data. One very widely utilized technique for 
training a neural network is a backpropagation training 
algorithm. However, there are other types of algorithms that 
can be utilized to set the "weights" in the network. 

[0015] When a large amount of steady-state data is avail- 
able to a network, the stored representation can be accurately 
modeled. However, some plants have a large amount of 
dynamic information associated therewith. This dynamic 
information reflects the fact that the inputs to the plant 
undergo a change which results in a corresponding change 
in the output. If a user desired to predict the final steady-state 
value of the plant, plant dynamics may not be important and 
this data could be ignored. However, there are situations 
wherein the dynamics of the plant are important during the 
prediction. It may be desirable to predict the path that an 
output will take from a beginning point to an end point. For 
example, if the input were to change in a step function from 
one value to another, a steady-state model that was accu- 
rately trained would predict the final steady-state value with 
some accuracy. However, the path between the starting point 
and the end point would not be predicted, as this would be 
subject to the dynamics of the plant, Further, in some control 
applications, it may be desirable to actually control the plant 
such that the plant dynamics were "constrained," this requir- 
ing some knowledge of the dynamic operation of the plant, 

[0016] In some applications, the actual historical data that 
is available as the training set has associated therewith a 
considerable amount of dynamic information. If the training 
data set had a large amount of steady-state information, an 
accurate model could easily be obtained for a steady-state 
model. However, if the historical data had a large amount of 
dynamic information associated therewith, i.e., the plant 
were not allowed to come to rest for a given data point, then 
there would be an error associated with the training opera- 
tion that would be a result of this dynamic component in the 
training data. This is typically the case for small data sets. 
This dynamic component must therefore be dealt with for 
small training data sets when attempting to train a steady- 
state model. 



[0017] When utilizing a model for the purpose of optimi- 
zation, it is necessary to train a model on one set of input 
values to predict another set of input values at future time. 
This will typically require a steady-state modeling tech- 
nique. In optimization, especially when used in conjunction 
with a control system, the optimization process will take a 
desired set of set points and optimizes those set points. 
However, these models are typically selected for accurate 
gain, a problem arises whenever the actual plant changes due 
to external influences, such as outside temperature, build up 
of slag, etc. Of course, one could regenerate the model with 
new parameters. However, the typical method is to actually 
measure the output of the plant, compare it with a predicted 
value to generate a "biased" value which sets forth the error 
in the plant as opposed to the model. This error is then 
utilized to bias the optimization network. However, to date 
this technique has required the use of steady-state models 
which are generally off-line models. The reason for this is 
that the actual values must "settle out" to reach a steady-state 
value before the actual bias can be determined. During 
operation of a plant, the outputs are dynamic. 

SUMMARY OF THE INVENTION 

[0018] The present invention disclosed herein comprises a 
method for modifying the values of select ones of a plurality 
of inputs to a plant. The plurality of inputs is received and 
then only N select ones of the received inputs are processed. 
These N select inputs are mapped in a first step of mapping 
through a first relationship. This first relationship defines P 
intermediate inputs numbering less than N, with the first 
relationship defining the relationship between the N select 
received inputs and the P intermediate inputs as a set of P 
intermediate values. These P intermediate values are then 
processed through a modifying process which modifies the 
P intermediate values in accordance with a predetermined 
modification algorithm to provide P modified intermediate 
inputs with corresponding P modified intermediate values. A 
second step of mapping is then performed on the P modified 
intermediate inputs and associated P modified intermediate 
values through the inverse of the first relationship. This 
provides N modified select input values at the output of the 
second mapping step that correspond to the N select 
received input values. 

[0019] In another aspect of the present invention, there is 
provided a post processing step on select ones of the N 
modified select input values. Each of the processed ones of 
the N modified select input values is processed with a 
predetermined process algorithm. Each of these process 
algorithms can be unique to the associated one of the N 
modified select input. 

BRIEF DESCRIPTION OF THE DRAWINGS 

[0020] For a more complete understanding of the present 
invention and the advantages thereof, reference is now made 
to the following description taken in conjunction with the 
accompanying Drawings in which: 

[0021] FIG. 1 illustrates a prior art Hammerstein model; 

[0022] FIG. 2 illustrates a block diagram of a modeling 
technique utilizing steady-state gain to define the gain of the 
dynamic model; 

[0023] FIGS. 3a-3d illustrate timing diagrams for the 
various outputs of the system of FIG. 2; 
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[0024] FIG. 4 illustrates a detailed block diagram of a 
dynamic model; 

[0025] FIG. 5 illustrates a block diagram of the operation 
of the. model of FIG. 4; 

[0026] FIG. 6 illustrates an example of the modeling 
technique of the embodiment of FIG. 2 utilized in a control 
environment; 

[0027] FIG. 7 illustrates a diagrammatic view of a change 
between two steady-state values; 

[0028] FIG. 8 illustrates a diagrammatic view of the 
approximation algorithm for changes in the steady-state 
value; 

[0029] FIG. 9 illustrates a block diagram of the dynamic 
model; 

[0030] FIG. 10 illustrates a detail of the control network 
utilizing an error constraining algorithm; 

[0031] FIGS, lla and lib illustrate plots of the input and 
output during optimization; 

[0032] FIG. 12 illustrates a plot depicting desired and 
predicted behavior; 

[0033] FIG. 13 illustrates various plots for controlling a 
system to force the predicted behavior to the desired behav- 
ior; 

[0034] FIG. 14 illustrates a plot of a trajectory weighting 
algorithm; 

[0035] FIG. 15 illustrates a plot for one type of constrain- 
ing algorithm; 

[0036] FIG. 16 illustrates a plot of the error algorithm as 
a function of time; 

[0037] FIG. 17 illustrates a flowchart depicting the statis- 
tical method for generating the filter and defining the end 
point for the constraining algorithm of FIG. 15; 

[0038] FIG. 18 illustrates a diagrammatic view of the 
optimization process; 

[0039] FIG. 18a illustrates a diagrammatic representation 
of the manner in which the path between steady-state values 
is mapped through the input and output space; 

[0040] FIG. 19 illustrates a flowchart for the optimization 
procedure; 

[0041] FIG. 20 illustrates a diagrammatic view of the 
input space and the error associated therewith; 

[0042] FIG. 21 illustrates a diagrammatic view of the 
confidence factor in the input space; 

[0043] FIG. 22 illustrates a block diagram of the method 
for utilizing a combination of a non- linear system and a first 
principals system; 

[0044] FIG. 23 illustrates an alternate embodiment of the 
embodiment of FIG. 22; 

[0045] FIG. 24 illustrates a plot of a pair of data with a 
defined delay associated therewith; 

[0046] FIG. 25 illustrates a diagrammatic view of the 
binning method for determining statistical independence; 



[0047] FIG. 26 illustrates a block diagram of a training 
method wherein delay is determined by statistical analysis; 

[0048] FIG. 27 illustrates a flow chart of the method for 
determining delays based upon statistical methods; 

[0049] FIG. 28 illustrates a prior art Weiner model; 

[0050] FIG. 29 illustrates a block diagram of a training 
method utilizing the system dynamics; 

[0051] FIG. 30 illustrates plots of input data, actual output 
data, and the filtered input data which has the plant dynamics 
impressed thereupon; 

[0052] FIG. 31 illustrates a flow chart for the training 
operation; 

[0053] FIG. 32 illustrates a diagrammatic view of the step 
test; 

[0054] FIG. 33 illustrates a diagrammatic view of a single 
step for u(t) and u(t); 

[0055] FIG. 34 illustrates a diagrammatic view of the 
pre-filter operation during training; 

[0056] FIG. 35 illustrates a diagrammatic view of a 
MIMO implementation of the training method of the present 
invention; and 

[0057] FIG. 36 illustrates a non-fully connected network; 

[0058] FIG. 37 illustrates a graphical user interface for 
selecting ranges of values for the dynamic inputs in order to 
train the dynamic model; 

[0059] FIG. 38 illustrates a flowchart depicting the selec- 
tion of data and training of the model; 

[0060] FIGS. 39 and 40 illustrate graphical user inter- 
faces for depicting both the actual historical response and 
the predictive response; 

[0061] FIG. 41 illustrates a block diagram of a predictive 
control system with a GUI interface; 

[0062] FIGS. 41-45 illustrate screen views for changing 
the number of variables that can be displayed from a given 
set; 

[0063] FIG. 46 illustrates a diagrammatic view of a plant 
utilizing on-line optimization; 

[0064] FIG. 47 illustrates a block diagram of the opti- 
mizer; 

[0065] FIG. 48 illustrates a plot of manipulatable vari- 
ables and controlled variables or outputs; 

[0066] FIGS. 49-51 illustrate plots of the dynamic opera- 
tion of the system and the bias; 

[0067] FIG. 52 illustrates a block diagram of a prior art 
optimizer utilizing steady-state; 

[0068] FIG. 53 illustrates a diagrammatic view for deter- 
mining the computed disturbance variables; 

[0069] FIG. 53a illustrates a block diagram of a residual 
activation neural network; 

[0070] FIG. 54 illustrates a block diagram for a steady 
state model utilizing the computer disturbance variables; 
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[0071] FIG, 55 illustrates an overall block diagram of an 
optimization circuit utilizing computed disturbance vari- 
ables; 

[0072] FIG. 56 illustrates a diagrammatic view of furnace/ 
boiler system which has associated therewith multiple levels 
of coal firing; 

[0073] FIG. 57 illustrates a top sectional view of the 
tangentially fired furnace; 

[0074] FIG. 57a illustrates a side cross-sectional view of 
the tangentially fired furnace; 

[0075] FIG. 58 illustrates a block diagram of one appli- 
cation of the on-line optimizer; 

[0076] FIG. 59 illustrates a block diagram of training 
algorithm for training a model using a multiple to single MV 
algorithm; and 

[0077] FIGS. 60 and 61 illustrate more detailed block 
diagrams of the embodiment of FIG. 57. 

[0078] FIG. 62 illustrates a conceptual diagrammatic view 
of the process for converting multiple inputs to a single input 
for optimization thereof; and 

[0079] FIG. 63 illustrates a more detailed block diagram 
of the structure in FIG. 62. 

DETAILED DESCRIPTION OF THE 
INVENTION 

[0080] Referring now to FIG. 1, there is illustrated a 
diagrammatic view of a Hammerstein model of the prior art. 
This is comprised of a non-linear static operator model 10 
and a linear dynamic model 12, both disposed in a series 
configuration. The operation of this model is described in H. 
T Su, and T J. McAvoy, "Integration of Multilayer Percep- 
tron Networks and Linear Dynamic Models: A Hammerstein 
Modeling Approach" to appear in I & EC Fundamentals, 
paper dated Jul. 7, 1992, which reference is incorporated 
herein by reference. Hammerstein models in general have 
been utilized in modeling non-linear systems for some time. 
The structure of the Hammerstein model illustrated in FIG. 
1 utilizes the non-linear static operator model 10 to trans- 
form the input U into intermediate variables H. The non- 
linear operator is usually represented by a finite polynomial 
expansion. However, this could utilize a neural network or 
any type of compatible modeling system. The linear 
dynamic operator model 12 could utilize a discreet dynamic 
transfer function representing the dynamic relationship 
between the intermediate variable H and the output Y. For 
multiple input systems, the non-linear operator could utilize 
a multilayer neural network, whereas the linear operator 
could utilize a two layer neural network. A neural network 
for the static operator is generally well known and described 
in U.S. Pat. No. 5,353,207, issued Oct. 4, 1994, and assigned 
to the present assignee, which is incorporated herein by 
reference. These type of networks are typically referred to as 
a multilayer feed-forward network which utilizes training in 
the form of back-propagation. This is typically performed on 
a large set of training data. Once trained, the network has 
weights associated therewith, which are stored in a separate 
database. 

[0081] Once the steady-state model is obtained, one can 
then choose the output vector from the hidden layer in the 



neural network as the intermediate variable for the Ham- 
merstein model. In order to determine the input for the linear 
dynamic operator, u(t), it is necessary to scale the output 
vector h(d) from the non-linear static operator model 10 for 
the mapping of the intermediate variable h(t) to the output 
variable of the dynamic model y(t), which is determined by 
the linear dynamic model. 

[0082] During the development of a linear dynamic model 
to represent the linear dynamic operator, in the Hammerstein 
model, it is important that the steady-state non-linearity 
remain the same. To achieve this goal, one must train the 
dynamic model subject to a constraint so that the non- 
linearity learned by the steady-state model remains 
unchanged after the training. This results in a dependency of 
the two models on each other. 

[0083] Referring now to FIG. 2, there is illustrated a block 
diagram of the modeling method in one embodiment, which 
is referred to as a systematic modeling technique. The 
general concept of the systematic modeling technique in the 
present embodiment results from the observation that, while 
process gains (steady-state behavior) vary with U's and 
Y's,(i.e., the gains are non-linear), the process dynamics 
seemingly vary with time only, (i.e., they can be modeled as 
locally linear, but time-varied). By utilizing non-linear mod- 
els for the steady-state behavior and linear models for the 
dynamic behavior, several practical advantages result. They 
are as follows: 

[0084] 1 . Completely rigorous models can be utilized 
for the steady-state part. This provides a credible 
basis for economic optimization. 

[0085] 2. The linear models for the dynamic part can 
be updated on-line, i.e., the dynamic parameters that 
are known to be time-varying can be adapted slowly. 

[0086] 3. The gains of the dynamic models and the 
gains of the steady-state models can be forced to be 
consistent (k=K). 

[0087] With further reference to FIG. 2, there are provided 
a static or steady-state model 20 and a dynamic model 22. 
The static model 20, as described above, is a rigorous model 
that is trained on a large set of steady-state data. The static 
model 20 will receive a procesa jnput U and provid e a 
predicted output Y. These are essentially steady-state values. 
The steady-state values at a given time are latched in various 
latches, an input latch 24 and an output latch 26. The latch 
24 contains the steady-state value of the input U„, and the 
latch 26 contains the steady-state output value Y M . The 
dynamic model 22 is utilized to predict the behavior of the 
plant when a change is made from a steady-state value of Y ss 
to a new value Y. The dynamic model 22 receives on the 
input the dynamic input value u and outputs a predicted 
dynamic value y. The value u is comprised of the difference 
between the new value U and the steady-state value in the 
latch 24, U sg . This is derived from a subtraction circuit 30 
which receives on the positive input thereof the output of the 
latch 24 and on the negative input thereof the new value of 
U. This therefore represents the delta change from the 
steady-state. Similarly, on the output the predicted overall 
dynamic value will be the sum of the output value of the 
dynamic model, y, and the steady-state output value stored 
in the latch 26, Y M These two values are summed with a 
summing block 34 to provide a predicted output Y. The 
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difference between the value output by the summing junc- 
tion 34 and the predicted value output by the static model 20 
is that the predicted value output by the summing junction 
20 accounts for the dynamic operation of the system during 
a change. For example, to process the input values that are 
in the input vector U by the static model 20, the rigorous 
model, can take significantly more time than running a 
relatively simple dynamic model. The method utilized in the 
present embodiment is to force the gain of the dynamic 
model 22 k d to equal the gain of the static model 20. 

[0088] In the static model 20, there is provided a storage 
block 36 which contains the static coefficients associated 
with the static model 20 and also the associated gain value 
K^. Similarly, the dynamic model 22 has a storage area 38 
that is operable to contain the dynamic coefficients and the 
gain value k d . One of the important aspects of the present 
embodiment is a link block 40 that is operable to modify the 
coefficients in the storage area 38 to force the value of k d to 
be equal to the value of K^. Additionally, there is an 
approximation block 41 that allows approximation of the 
dynamic gain k d between the modification updates. 

[0089] Systematic Model 

[0090] The linear dynamic model 22 can generally be 
represented by the following equations: 



n n 

6y(t) = Yu ~ d - " Z a ' t5y{t ~ 0 



(007) 



[0091] where: 
&«(?)="(')-«« 



9 



(008) 
(009) 

[0092] and t is time, aj and b { are real numbers, d is a time 
delay, u(t) is an input and y(t) an 



"(B) 



= k = 



(10) 



1 + 1 fl/tf- 1 



[0093] where B is the backward shift operator B(x(t))» 
x(t-l), t=time, the a { and b { are real numbers, I is the number 
of discreet time intervals in the dead-time of the process, and 
n is the order of the model. This is a general representation 
of a linear dynamic model, as contained in George E. P. Box 
and G. M. Jenkins, "TIME SERIES ANALYSIS forecasting 
and control", Holden-Day, San Francisco, 1976, Section 
10.2, Page 345. This reference is incorporated herein by 
reference. 

[0094] The gain of this model can be calculated by setting 
the value of B equal to a value of "1". The gain will then be 
defined by the following equation: 



mL =kd= 



1 + £ Qi 



(11) 



[0095] The a t contain the dynamic signature of the pro- 
cess, its unforced, natural response characteristic. They are 
independent of the process gain. The b. contain part of the 
dynamic signature of the process; however, they alone 
contain the result of the forced response. The bj determine 
the gain k of the dynamic model. See: J. L. Shearer, A. T. 
Murphy, and H. H. Richardson, "Introduction to System 
Dynamics", Addison-Wesley, Reading, Massachusetts, 
1967, Chapter 12. This reference is incorporated herein by 
reference. 

[0096] Since the gain of the steady-state model is 
known, the gain k d of the dynamic model can be forced to 
match the gain of the steady-state model by scaling the 
parameters. The values of the static and dynamic gains are 
set equal with the value of b. scaled by the ratio of the two 
gains: 



(12) 
(13) 



[0097] This makes the dynamic model consistent with its 
steady-state counterpart. Therefore, each time the steady- 
state value changes, this corresponds to a gain K„ of the 
steady-state model. This value can then be utilized to update 
the gain k d of the dynamic model and, therefore, compensate 
for the errors associated with the dynamic model wherein 
the value of k d is determined based on perturbations in the 
plant on a given set of operating conditions. Since all 
operating conditions are not modeled, the step of varying the 
gain will account for changes in the steady-state starting 
points. 

[0098] Referring now to FIGS. 3a-3</, there are illustrated 
plots of the system operating in response to a step function 
wherein the input value U changes from a value of 100 to a 
value of 110, In FIG. 3a, the value of 100 is referred to as 
the previous steady-state value U aa . In FIG. 36, the value of 
u varies from a value of 0 to a value of 10, this representing 
the delta between the steady-state value of U Bg to the level of 
110, represented by reference numeral 42 in FIG. 3a. 
Therefore, in FIG. 36 the value of u will go from 0 at a level 
44, to a value of 10 at a level 46. In FIG. 3c, the output Y 
is represented as having a steady-state value Y„ of 4 at a 
level 48. When the input value U rises to the level 42 with 
a value of 110, the output value will rise. This is a predicted 
value. The predicted value which is the proper output value 
is represented by a level 50, which level 50 is at a value of 
5. Since the steady-state value is at a value of 4, this means 
that the dynamic system must predict a difference of a value 
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of 1. This is represented by FIG. 3d wherein the dynamic 
output value y varies from a level 54 having a value of 0 to 
a level 56 having a value of 1.0. However, without the gain 
scaling, the dynamic model could, by way of example, 
predict a value for y of 1.5, represented by dashed level 58, 
if the steady-state values were outside of the range in which 
the dynamic model was trained. This would correspond to a 
value of 5.5 at a level 60 in the plot of FIG. 3c. It can be seen 
that the dynamic model merely predicts the behavior of the 
plant from a starting point to a stopping point, not taking into 
consideration the steady-state values. It assumes that the 
steady-state values are those that it was trained upon. If the 
gain k d were not scaled, then the dynamic model would 
assume that the steady-state values at the starting point were 
the same that it was trained upon. However, the gain scaling 
link between the steady-state model and the dynamic model 
allow the gain to be scaled and the parameter bj to be scaled 
such that the dynamic operation is scaled and a more 
accurate prediction is made which accounts for the dynamic 
properties of the system. 

[0099] Referring now to FIG. 4, there is illustrated a block 
diagram of a method for determining the parameters a^ b £ . 
This is usually achieved through the use of an identification 
algorithm, which is conventional. This utilizes the (u(t),y(t)) 
pairs to obtain the aj and b t parameters. In the preferred 
embodiment, a recursive identification method is utilized 
where the a { and parameters are updated with each new 
( u i(0> v i(0) P a * r See: T. Eykhoff, "System Identification", 
John Wiley & Sons, New York, 1974, Pages 38 and 39, et. 
seq., and H. Kurz and W. Godecke, "Digital Parameter- 
Adaptive Control Processes with Unknown Dead Time", 
Automatica, Vol. 17, No. 1, 1981, pp. 245-252, which 
references are incorporated herein by reference. 

[0100] In the technique of FIG. 4, the dynamic model 22 
has the output thereof input to a parameter-adaptive control 
algorithm block 60 which adjusts the parameters in the 
coefficient storage block 38, which also receives the scaled 
values of k, fy. This is a system that is updated on a periodic 
basis, as defined by timing block 62. The control algorithm 
60 utilizes both the input u and the output y for the purpose 
of determining and updating the parameters in the storage 
area 38. 

[0101] Referring now to FIG. 5, there is illustrated a block 
diagram of the preferred method. The program is initiated in 
a block 68 and then proceeds to a function block 70 to update 
the parameters a £ , b { utilizing the (u(I),y(I)) pairs. Once these 
are updated, the program flows to a function block 72 
wherein the steady-state gain factor K is received, and then 
to a function block 74 to set the dynamic gain to the 
steady-state gain, i.e., provide the scaling function described 
hereinabove. This is performed after the update. This pro- 
cedure can be used for on-line identification, non-linear 
dynamic model prediction and adaptive control. 

[0102] Referring now to FIG. 6, there is illustrated a block 
diagram of one application of the present embodiment 
utilizing a control environment. A plant 78 is provided which 
receives input values u(t) and outputs an output vector y(t). 
The plant 78 also has measurable state variables s(t). A 
predictive model 80 is provided which receives the input 
values u(t) and the slate variables s(t) in addition to the 
output value y(t). The steady-state model 80 is operable to 
output a predicted value of both y(t) and also of a future 



input value u(t+l). This constitutes a steady-state portion of 
the system. The predicted steady-state input value is U 8S with 
the predicted steady-state output value being Y ss . In a 
conventional control scenario, the steady-state model 80 
would receive as an external input a desired value of the 
output y d (t) which is the desired value that the overall 
control system seeks to achieve. This is achieved by con- 
trolling a distributed control system (DCS) 86 to produce a 
desired input to the plant. This is referred to as u(t+l), a 
future value. Without considering the dynamic response, the 
predictive model 80, a steady-state model, will provide the 
steady-state values. However, when a change is desired, this 
change will effectively be viewed as a "step response". 

[0103] To facilitate the dynamic control aspect, a dynamic 
controller 82 is provided which is operable to receive the 
input u(t), the output value y(t) and also the steady-state 
values U„ and Y ss and generate the output u(t+l). The 
dynamic controller effectively generates the dynamic 
response between the changes, i.e., when the steady-state 
value changes from an initial steady-state value U sa l , Y 1 ^ to 
a final steady-state value U f ss) Y f ss . 

[0104] During the operation of the system, the dynamic 
controller 82 is operable in accordance with the embodiment 
of FIG. 2 to update the dynamic parameters of the dynamic 
controller 82 in a block 88 with a gain link block 90, which 
utilizes the value K„ from a steady-state parameter block in 
order to scale the parameters utilized by the dynamic con- 
troller 82, again in accordance with the above described 
method. In this manner, the control function can be realized. 
In addition, the dynamic controller 82 has the operation 
thereof optimized such that the path traveled between the 
initial and final steady-state values is achieved with the use 
of the optimizer 83 in view of optimizer constraints in a 
block 85. In general, the predicted model (steady-state 
model) 80 provides a control network function that is 
operable to predict the future input values. Without the 
dynamic controller 82, this is a conventional control net- 
work which is generally described in U.S. Pat. No. 5,353, 
207, issued Oct. 4, 1994, to the present assignee, which 
patent is incorporated herein by reference. 

[0105] Approximate Systematic Modeling 

[0106] For the modeling techniques described thus far, 
consistency between the steady-state and dynamic models is 
maintained by rescaling the b { parameters at each time step 
utilizing equation 13. If the systematic model is to be 
utilized in a Model Predictive Control (MPC) algorithm, 
maintaining consistency may be computationally expensive. 
These types of algorithms are described in C. E. Garcia, D. 
M. Prett and M. Morari. Model predictive control: theory 
and practice — a survey, Automatica, 25:335-348,1989; D. E. 
Seborg, T. F. Edgar, and D. A. Mellichamp. Process Dynam- 
ics and Control. John Wiley and Sons, New York, N.Y, 
1989. These references are incorporated herein by reference. 
For example, if the dynamic gain kj is computed from a 
neural network steady-state model, it would be necessary to 
execute the neural network module each time the model was 
iterated in the MPC algorithm. Due to the potentially large 
number of model iterations for certain MPC problems, it 
could be computationally expensive to maintain a consistent 
model. In this case, it would be better to use an approximate 
model which does not rely on enforcing consistencies at 
each iteration of the model. 



06/03/2004, EAST version: 1.4.1 



US 2003/0018399 Al 



Jan. 23, 2003 



[0107] Referring now to FIG. 7, there is illustrated a 
diagram for a change between steady-state values. As illus- 
trated, the steady-state model will make a change from a 
steady-state value at a line 100 to a steady-state value at a 
line 102. A transition between the two steady-state values 
can result in unknown settings. The only way to insure that 
the settings for the dynamic model between the two steady- 
state values, an initial steady-state value K^ 1 and a final 
steady-state gain K,/, would be to utilize a step operation, 
wherein the dynamic gain k d was adjusted at multiple 
positions during the change. However, this may be compu- 
tationally expensive. As will be described hereinbelow, an 
approximation algorithm is utilized for approximating the 
dynamic behavior between the two steady-state values uti- 
lizing a quadratic relationship. This is defined as a behavior 
line 104, which is disposed between an envelope 106, which 
behavior line 104 will be described hereinbelow. 

[0108] Referring now to FIG. 8, there is illustrated a 
diagrammatic view of the system undergoing numerous 
changes in steady-state value as represented by a stepped 
line 108. The stepped line 108 is seen to vary from a first 
steady-state value at a level 110 to a value at a level 112 and 
then down to a value at a level 114, up to a value at a level 
116 and then down to a final value at a level 118. Each of 
these transitions can result in unknown states. With the 
approximation algorithm that will be described hereinbelow, 
it can be seen that, when a transition is made from level 110 
to level 112, an approximation curve for the dynamic 
behavior 120 is provided. When making a transition from 
level 114 to level 116, an approximation gain curve 124 is 
provided to approximate the steady-state gains between the 
two levels 114 and 116. For making the transition from level 
116 to level 118, an approximation gain curve 126 for the 
steady-state gain is provided. It can therefore be seen that the 
approximation curves 120-126 account for transitions 
between steady-state values that are determined by the 
network, it being noted that these are approximations which 
primarily maintain the steady-state gain within some type of 
error envelope, the envelope 106 in FIG. 7. 

[0109] The approximation is provided by the block 41 
noted in FIG. 2 and can be designed upon a number of 
criteria, depending upon the problem that it will be utilized 
to solve. The system in the preferred embodiment, which is 
only one example, is designed to satisfy the following 
criteria: 

[0110] 1. Computational Complexity: The approxi- 
mate systematic model will be used in a Model 
Predictive Control algorithm, therefore, it is required 
to have low computational complexity. 

[0111] 2. Localized Accuracy: The steady-state 
model is accurate in localized regions. These regions 
represent the steady-state operating regimes of the 
process. The steady-state model is significantly less 
accurate outside these localized regions. 

[0112] 3. Final Steady-State: Given a steady-state set 
point change, an optimization algorithm which uses 
the steady-state model will be used to compute the 
steady-state inputs required to achieve the set point. 
Because of item 2, it is assumed that the initial and 
final steady-states associated with a set-point change 
are located in regions accurately modeled by the 
steady-state model. 



[0113] Given the noted criteria, an approximate systematic 
model can be constructed by enforcing consistency of the 
steady-state and dynamic model at the initial and final 
steady-state associated with a set point change and utilizing 
a linear approximation at points in between the two steady- 
states. This approximation guarantees that the model is 
accurate in regions where the steady-state model is well 
known and utilizes a linear approximation in regions where 
the steady-state model is known to be less accurate. In 
addition, the resulting model has low computational com- 
plexity. For purposes of this proof, Equation 13 is modified 
as follows: 



biK^uU-d-l))] 



('•§«] 



(14) 



;=i 



[0114] This new equation 14 utilizes K as (u(t-d— 1)) instead 
of K^uO)) as the consistent gain, resulting in a systematic 
model which is delay invariant. 

[0115] The approximate systematic model is based upon 
utilizing the gains associated with the initial and final 
steady-slate values of a set-point change. The initial steady- 
state gain is denoted K'„ while the initial steady-state input 
is given by U 1 . The final steady-state gain is K f ss and the 
final input is U Sfl . Given these values, a linear approximation 
to the gain is given by: 



(15) 



[0116] Substituting this approximation into Equation 13 
and replacing u^-d-l)-^ by 6«(t-d-l) yields: 



XL) 



(16) 



-6u(t-d-i), 



hjjecltd = — + 2 / n "\ , " 



[0117] To simplify the expression, define the variable 
bj-Bar as: 



(17) 



M 

and gj as : 

^M + >>I(*£-KL) 



8j"- 



(18) 
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[0118] Equation 16 may be written as: 

6j,.«ied=fii+*M'-<*-0- (19) 

[0119] Finally, substituting the scaled b's back into the 
original difference Equation 7, the following expression for 
the approximate systematic model is obtained: 

n (20) 
6y(t) = 5|tf«(f-rf-0 + 
i-l 

n n 

£ g-M* -d- i 2 )5u(t - d - i) - £ a -M" ~ 0 



[0120] The linear approximation for gain results in a 
quadratic difference equation for the output. Given Equation 
20, the approximate systematic model is shown to be of low 
computational complexity. It may be used in a MPC algo- 
rithm to efficiently compute the required control moves for 
a transition from one steady-state to another after a set-point 
change. Note that this applies to the dynamic gain variations 
between steady-state transitions and not to the actual path 
values. 

[0121] Control System Error Constraints 

[0122] Referring now to FIG. 9, there is illustrated a block 
diagram of the prediction engine for the dynamic controller 
82 of FIG. 6. The prediction engine is operable to essentially 
predict a value of y(t) as the predicted future value y(t+l). 
Since the prediction engine must determine what the value 
of the output y(t) is at each future value between two 
steady-state values, it is necessary to perform these in a 
"step" manner. Therefore, there will be k steps from a value 
of zero to a value of N, which value at k-N is the value at 
the "horizon", the desired value. This, as will be described 
hereinbelow, is an iterative process, it being noted that the 
terminology for "(t+1)" refers to an incremental step, with 
an incremental step for the dynamic controller being smaller 
than an incremented step for the steady-state model. For the 
steady-state model, "y(t+N)" for the dynamic model will be, 
"y(t+l)" for the steady state The value y(t+l) is defined as 
follows: 

y(t+*l)=ay{i)+aj,y(t-\)+biu{t-d-l}+b?. u(t-d-2) (021) 

[0123] With further reference to FIG. 9, the input values 
u(t) for each (u,y) pair are input to a delay line 140. The 
output of the delay line provides the input value u(t) delayed 
by a delay value "d". There are provided only two operations 
for multiplication with the coefficients b A and b 2 , such that 
only two values u(t) and u(t-l) are required. These are both 
delayed and then multiplied by the coefficients b 1 and b 2 and 
then input to a summing block 141. Similarly, the output 
value y p (t) is input to a delay line 142, there being two 
values required for multiplication with the coefficients 3j 
and a 2 . The output of this multiplication is then input to the 
summing block 141. The input to the delay line 142 is either 
the actual input value y a (t) or the iterated output value of the 
summation block 141, which is the previous value computed 
by the dynamic controller 82. Therefore, the summing block 
141 will output the predicted value y(t+l) which will then be 
input to a multiplexor 144. The multiplexor 144 is operable 
to select the actual output y a (t) on the first operation and, 
thereafter, select the output of the summing block 141. 



Therefore, for a step value of k«0 the value y 8 (t) will be 
selected by the multiplexor 144 and will be latched in a latch 
145. The latch 145 will provide the predicted value y p (t+k) 
on an output 146. This is the predicted value of y(t) for a 
given k that is input back to the input of delay line 142 for 
multiplication with the coefficients ^ 1 and a^. This is iterated 
for each value of k from k=0 to k«N. 

[0124] The a a and values are fixed, as described above, 
with the bj and b 2 values scaled. This scaling operation is 
performed by the coefficient modification block 38. How- 
ever, this only defines the beginning steady-state value and 
the final steady-state value, with the dynamic controller and 
the optimization routines described in the present applica- 
tion defining how the dynamic controller operates between 
the steady -state values and also what the gain of the dynamic 
controller is. The gain specifically is what determines the 
modification operation performed by the coefficient modi- 
fication block 38. 

[0125] In FIG. 9, the coefficients in the coefficient modi- 
fication block 38 are modified as described hereinabove with 
the information that is derived from the steady-state model. 
The steady-state model is operated in a control application, 
and is comprised in part of a forward steady -state model 141 
which is operable to receive the steady-state input value 
U M (t) and predict the steady-state output value Y 8S (t) This 
predicted value is utilized in an inverse steady-state model 
143 to receive the desired value y d (t) and the predicted 
output of the steady-state model 141 and predict a future 
steady-state input value or manipulated value U M (t+N) and 
also a future steady-state input value Y 8S (t+N) in addition to 
providing the steady-state gain K^, As described herein- 
above, these are utilized to generate scaled b-values. These 
b -values are utilized to define the gain kj of the dynamic 
model. In can therefore be seen that this essentially takes a 
linear dynamic model with a fixed gain and allows it to have 
a gain thereof modified by a non-linear model as the 
operating point is moved through the output space. 

[0126] Referring now to FIG. 10, there is illustrated a 
block diagram of the dynamic controller and optimizer. The 
dynamic controller includes a dynamic model 149 which 
basically defines the predicted value y p (k) as a function of 
the inputs y(t), s(t) and u(t). This was essentially the same 
model that was described hereinabove with reference to 
FIG. 9. The model 149 predicts the output values y p (k) 
between the two steady-state values, as will be described 
hereinbelow. The model 149 is predefined and utilizes an 
identification algorithm to identify the a a , a 2 , b a and b 2 
coefficients during training. Once these are identified in a 
training and identification procedure, these are "fixed". 
However, as described hereinabove, the gain of the dynamic 
model is modified by scaling the coefficients b a and b 2 . This 
gain scaling is not described with respect to the optimization 
operation of FIG. 10 although it can be incorporated in the 
optimization operation. 

[0127] The output of model 149 is input to the negative 
input of a summing block 150. Summing block 150 sums the 
predicted output y p (k) with the desired output y^t). In effect, 
the desired value of y d (t) is effectively the desired steady- 
state value Y f as , although it can be any desired value. The 
output of the summing block 150 comprises an error value 
which is essentially the difference between the desired value 
y d (t) and the predicted value y p (k). The error value is 
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modified by an error modification block 151, as will be 
described hereinbelow, in accordance with error modifica- 
tion parameters in a block 152. The modified error value is 
then input to an inverse model 153, which basically per- 
forms an optimization routine to predict a change in the 
input value u(t). In effect, the optimizer 153 is utilized in 
conjunction with the model 149 to minimize the error output 
by summing block 150. Any optimization function can be 
utilized, such as a Monte Carlo procedure. However, in the 
present embodiment, a gradient calculation is utilized. In the 
gradient method, the gradient d(y)/d(u) is calculated and 
then a gradient solution performed as follows: 



[0128] The optimization function is performed by the 
inverse model 153 in accordance with optimization con- 
straints in a block 154. An iteration procedure is performed 
with an iterate block 155 which is operable to perform an 
iteration with the combination of the inverse model 153 and 
the predictive model 149 and output on an output line 156 
the future value u(t+k+l). For k=0, this will be the initial 
steady-state value and for k=N, this will be the value at the 
horizon, or at the next steady-state value. During the itera- 
tion procedure, the previous value of u(t+k) has the change 
value Au added thereto. This value is utilized for that value 
of k until the error is within the appropriate levels. Once it 
is at the appropriate level, the next u(t+k) is input to the 
model 149 and the value thereof optimized with the iterate 
block 155. Once the iteration procedure is done, it is latched. 
As will be described hereinbelow, this is a combination of 
modifying the error such that the actual error output by the 
block 150 is not utilized by the optimizer 153 but, rather, a 
modified error is utilized. Alternatively, different optimiza- 
tion constraints can be utilized, which are generated by the 
block 154, these being described hereinbelow. 

[0129] Referring now to FIGS. 11a and lib, there are 
illustrated plots of the output y(t+k) and the input u k (t+k+l), 
for each k from the initial steady-state value to the horizon 
steady-state value at k-N. With specific reference to FIG. 
11a, it can be seen that the optimization procedure is 
performed utilizing multiple passes. 

[0130] In the first pass, the actual value u fl (t+k) for each k 
is utilized to determine the values of y(t+k) for each u,y pair. 
This is then accumulated and the values processed through 
the inverse model 153 and the iterate block 155 to minimize 
the error. This generates a new set of inputs u k (t+k+l) 
illustrated in FIG. lib. Therefore, the optimization after 
pass 1 generates the values of u(t+k+l) for the second pass. 
In the second pass, the values are again optimized in 
accordance with the various constraints to again generate 
another set of values for u(t+k+l). This continues until the 
overall objective function is reached. This objective function 
is a combination of the operations as a function of the error 
and the operations as a function of the constraints, wherein 
the optimization constraints may control the overall opera- 
tion of the inverse model 153 or the error modification 
parameters in block 152 may control the overall operation. 
Each of the optimization constraints will be described in 
more detail hereinbelow. 



[0131] Referring now to FIG. 12, there is illustrated a plot 
of y d (t) and y p (t). The predicted value is represented by a 
waveform 170 and the desired output is represented by a 
waveform 172, both plotted over the horizon between an 
initial steady-state value Y l 6S and a final steady-state value 
. It can be seen that the desired waveform prior to k=0 
is substantially equal to the predicted output. At k=0, the 
desired output waveform 172 raises its level, thus creating 
an error. It can be seen that at k»0, the error is large and the 
system then must adjust the manipulated variables to mini- 
mize the error and force the predicted value to the desired 
value. The objective function for the calculation of error is 
of the form: 



j k 



[0132] where: 

[0133] Du u is the change in input variable (IV) I at 
time interval 1 

[0134] A, is the weight factor for control variable 
(CV)j 

[0135] y p (t) is the predicted value of CV j at time 
interval k 

[0136] y d (t) is the desired value of CV j. 
[0137] Trajectory Weighting 

[0138] The present system utilizes what is referred to as 
"trajectory weighting" which encompasses the concept that 
one does not put a constant degree of importance on the 
future predicted process behavior matching the desired 
behavior at every future time set, i.e., at low k- values. One 
approach could be that one is more tolerant of error in the 
near term (low k-values) than farther into the future (high 
k-values). The basis for this logic is that the final desired 
behavior is more important than the path taken to arrive at 
the desired behavior, otherwise the path traversed would be 
a step function. This is illustrated in FIG. 13 wherein three 
possible predicted behaviors are illustrated, one represented 
by a curve 174 which is acceptable, one is represented by a 
different curve 176, which is also acceptable and one rep- 
resented by a curve 178, which is unacceptable since it goes 
above the desired level on curve 172. Curves 174-178 define 
the desired behavior over the horizon for k-1 to N. 

[0139] In Equation 23, the predicted curves 174-178 
would be achieved by forcing the weighting factors Aj to be 
time varying. This is illustrated in FIG. 14. In FIG. 14, the 
weighting factor A as a function of time is shown to have an 
increasing value as time and the value of k increases. This 
results in the errors at the beginning of the horizon (low 
k-values) being weighted much less than the errors at the end 
of the horizon (high k-values). The result is more significant 
than merely redistributing the weights out to the end of the 
control horizon at k«N. This method also adds robustness, or 
the ability to handle a mismatch between the process and the 
prediction model. Since the largest error is usually experi- 
enced at the beginning of the horizon, the largest changes in 
the independent variables will also occur at this point. If 
there is a mismatch between the process and the prediction 
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(model error), these initial moves will be large and some- 
what incorrect, which can cause poor performance and 
eventually instability. By utilizing the trajectory weighting 
method, the errors at the beginning of the horizon are 
weighted less, resulting in smaller changes in the indepen- 
dent variables and, thus, more robustness. 

[0140] Error Constraints 

[0141] Referring now to FIG. 15, there are illustrated 
constraints that can be placed upon the error. There is 
illustrated a predicted curve 180 and a desired curve 182, 
desired curve 182 essentially being a flat line. It is desirable 
for the error between curve 180 and 182 to be minimized. 
Whenever a transient occurs at t — 0, changes of some sort 
will be required. It can be seen that prior to t=0, curve 182 
and 180 are substantially the same, there being very little 
error between the two. However, after some type of transi- 
tion, the error will increase. If a rigid solution were utilized, 
the system would immediately respond to this large error 
and attempt to reduce it in as short a time as possible. 
However, a constraint frustum boundary 184 is provided 
which allows the error to be large at t=0 and reduces it to a 
minimum level at a point 186. At point 186, this is the 
minimum error, which can be set to zero or to a non-zero 
value, corresponding to the noise level of the output variable 
to be controlled. This therefore encompasses the same 
concepts as the trajectory weighting method in that final 
future behavior is considered more important that near term 
behavior. The ever shrinking minimum and/or maximum 
bounds converge from a slack position at t~0 to the actual 
final desired behavior at a point 186 in the constraint frustum 
method. 

[0142] The difference between constraint frustum and 
trajectory weighting is that constraint frustums are an abso- 
lute limit (hard constraint) where any behavior satisfying the 
limit is just as acceptable as any other behavior that also 
satisfies the limit. Trajectory weighting is a method where 
differing behaviors have graduated importance in time. It 
can be seen that the constraints provided by the technique of 
FIG. 15 requires that the value y p (t) is prevented from 
exceeding the constraint value. Therefore, if the difference 
between y d (t) and y p (t) is greater than that defined by the 
constraint boundary, then the optimization routine will force 
the input values to a value that will result in the error being 
less than the constraint value. In effect, this is a "clamp" on 
the difference between y p (t) and y d (t). In the trajectory 
weighting method, there is no "clamp*' on the difference 
therebetween; rather, there is merely an attenuation factor 
placed on the error before input to the optimization network. 

[0143] Trajectory weighting can be compared with other 
methods, there being two methods that will be described 
herein, the dynamic matrix control (DMC) algorithm and the 
identification and command (IdCom) algorithm. The DMC 
algorithm utilizes an optimization to solve the control prob- 
lem by minimizing the objective function: 



srZ Z (^*(/w-/(o)+z *i*X (Af/ '-' )2 

j k i 1 



[0144] where B { is the move suppression factor for input 
variable I. This is described in Cutler, C. R. and B. L. 



Ramaker, Dynamic Matrix Control — A Computer Control 
Algorithm, AIChE National Meeting, Houston, Tex. (April, 
1979), which is incorporated herein by reference. 

[0145] It is noted that the weights Aj and desired values 
y d (t) are constant for each of the control variables. As can be 
seen from Equation 24, the optimization is a trade off 
between minimizing errors between the control variables 
and their desired values and minimizing the changes in the 
independent variables. Without the move suppression term, 
the independent variable changes resulting from the set point 
changes would be quite large due to the sudden and imme- 
diate error between the predicted and desired values. Move 
suppression limits the independent variable changes, but for 
all circumstances, not just the initial errors. 

[0146] The IdCom algorithm utilizes a different approach. 
Instead of a constant desired value, a path is defined for the 
control variables to take from the current value to the desired 
value. This is illustrated in FIG. 16. This path is a more 
gradual transition from one operation point to the next. 
Nevertheless, it is still a rigidly defined path that must be 
met. The objective function for this algorithm takes the 
form: 



[0147] This technique is described in Richalet, J. A. Rault, 
J. L. Testud, and J. Papon, Model Predictive Heuristic 
Control: Applications to Industrial Processes, Automatica, 
14, 413-428 (1978), which is incorporated herein by refer- 
ence. It should be noted that the requirement of Equation 25 
at each time interval is sometimes difficult. In fact, for 
control variables that behave similarly, this can result in 
quite erratic independent variable changes due to the control 
algorithm attempting to endlessly meet the desired path 
exactly. 

[0148] Control algorithms such as the DMC algorithm that 
utilize a form of matrix inversion in the control calculation, 
cannot handle control variable hard constraints directly. 
They must treat them separately, usually in the form of a 
steady-state linear program. Because this is done as a 
steady-state problem, the constraints are time invariant by 
definition. Moreover, since the constraints are not part of a 
control calculation, there is no protection against the con- 
troller violating the hard constraints in the transient while 
satisfying them at steady-state. 

[0149] With further reference to FIG. 15, the boundaries 
at the end of the envelope can be defined as described 
hereinbelow. One technique described in the prior art, W. 
Edwards Deming, "Out of the Crisis," Massachusetts Insti- 
tute of Technology, Center for Advanced Engineering Study, 
Cambridge Mass., Fifth Printing, September 1988, pages 
327-329, describes various Monte Carlo experiments that set 
forth the premise that any control actions taken to correct for 
common process variation actually may have a negative 
impact, which action may work to increase variability rather 
than the desired effect of reducing variation of the controlled 
processes. Given that any process has an inherent accuracy, 
there should be no basis to make a change based on a 
difference that lies within the accuracy limits of the system 



06/03/2004, EAST Version: 1.4.1 



US 2003/0018399 Al 



11 



Jan. 23, 2003 



utilized to control it. At present, commercial controllers fail 
to recognize the fact that changes are undesirable, and 
continually adjust the process, treating all deviation from 
target, no matter how small, as a special cause deserving of 
control actions, i.e., they respond to even minimal changes. 
Over adjustment of the manipulated variables therefore will 
result, and increase undesirable process variation. By plac- 
ing limits on the error with the present filtering algorithms 
described herein, only controller actions that are proven to 
be necessary are allowed, and thus, the process can settle 
into a reduced variation free from unmerited controller 
disturbances. The following discussion will deal with one 
technique for doing this, this being based on statistical 
parameters. 

[0150] Filters can be created that prevent model-based 
controllers from taking any action in the case where the 
difference between the controlled variable measurement and 
the desired target value are not significant. The significance 
level is defined by the accuracy of the model upon which the 
controller is statistically based. This accuracy is determined 
as a function of the standard deviation of the error and a 
predetermined confidence level. The confidence level is 
based upon the accuracy of the training. Since most training 
sets for a neural network-based model will have "holes" 
therein, this will result in inaccuracies within the mapped 
space. Since a neural network is an empirical model, it is 
only as accurate as the training data set. Even though the 
model may not have been trained upon a given set of inputs, 
it will extrapolate the output and predict a value given a set 
of inputs, even though these inputs are mapped across a 
space that is questionable. In these areas, the confidence 
level in the predicted output is relatively low. This is 
described in detail in U.S. patent application Ser. No. 
08/025,184, filed Mar. 2, 1993, which is incorporated herein 
by reference. 

[0151] Referring now to FIG. 17, there is illustrated a 
flowchart depicting the statistical method for generating the 
filter and defining the end point 186 in FIG. 15. The 
flowchart is initiated at a start block 200 and then proceeds 
to a function block 202, wherein the control values u(t+l) 
are calculated. However, prior to acquiring these control 
values, the filtering operation must be a processed. The 
program will flow to a function block 204 to determine the 
accuracy of the controller. This is done off-line by analyzing 
the model predicted values compared to the actual values, 
and calculating the standard deviation of the error in areas 
where the target is undisturbed. The model accuracy of e m (t) 
is defined as follows: 

^(0-«(0-pCO (026) 
[0152] where: 
[0153] e m -model error, 
[0154] a=actual value 
[0155] p=model predicted value 

[0156] The model accuracy is defined by the following 
equation: 

Acc-H+o m (027) 

[0157] where: 

[0158] Acc=accuracy in terms of minimal detector 
error 



H = significance level = 1 67% confidence 
= 2 95% confidence 
= 3 99.5% confidence 

[0159] a m =standard deviation of e m (t). 

[0160] The program then flows to a function block 206 to 
compare the controller error e c (t) with the model accuracy. 
This is done by taking the difference between the predicted 
value (measured value) and the desired value. This is the 
controller error calculation as follows: 

e c (t)-d(t)-m(t) (028) 

[0161] where: 

[0162] e c =controller error 

[0163] d«desired value 

[0164] m-measured value 

[0165] The program will then flow to a decision block 208 
to determine if the error is within the accuracy limits. The 
determination as to whether the error is within the accuracy 
limits is done utilizing Shewhart limits. With this type of 
limit and this type of filter, a determination is made as to 
whether the controller error e c (t) meets the following con- 
ditions: e c (t)^-l*Acc and e c (t)^+l*Acc, then either the 
control action is suppressed or not suppressed. If it is within 
the accuracy limits, then the control action is suppressed and 
the program flows along a "Y" path. If not, the program will 
flow along the "N" path to function block 210 to accept the 
u(t+l) values. If the error lies within the controller accuracy, 
then the program flows along the "Y" path from decision 
block 208 to a function block 212 to calculate the running 
accumulation of errors. This is formed utilizing a CUSUM 
approach. The controller CUSUM calculations are done as 
follows: 

5 low -min(0,5 low (r-l)+rf(0-/«(0)-2:(/«)+*) (029) 
S u -max(P, S u (/-lH^0-m(/))-2(j«)H) (030) 

[0166] where: 

[0167] S hl «Running Positive Qsum 

[0168] Si ow «Running Negative Qsum 

[0169] k=1\ining factor — minimal detectable 
change threshold 

[0170] with the following defined: 

[0171] Hq=significance level. Values of (j,k) can 
be found so that the CUSUM control chart will 
have significance levels equivalent to Shewhart 
control charts. 

[0172] The program will then flow to a decision block 214 
to determine if the CUSUM limits check out, i.e., it will 
determine if the Qsum values are within the limits. If the 
Qsum, the accumulated sum error, is within the established 
limits, the program will then flow along the "Y" path. And, 
if it is not within the limits, it will flow along the "N" path 
to accept the controller values u(t+l). The limits are deter- 
mined if both the value of S^+l'Hq and S low ^-l*Hq. 
Both of these actions will result in this program flowing 
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along the "Y" path. If it flows along the "N" path, the sum 
is set equal to zero and then the program flows to the 
function block 210. If the Qsum values are within the limits, 
it flows along the "Y" path to a function block 218 wherein 
a determination is made as to whether the user wishes to 
perturb the process. If so, the program will flow along the 
"Y" path to the function block 210 to accept the control 
values u(t +1), If not, the program will flow along the "N" 
path from decision block 218 to a function block 222 to 
suppress the controller values u(t +1). The decision block 
218, when it flows along the "Y" path, is a process that 
allows the user to re-identify the model for on-line adapta- 
tion, i.e., retrain the model. This is for the purpose of data 
collection and once the data has been collected, the system 
is then reactivated. 

[0173] Referring now to FIG. 18, there is illustrated a 
block diagram of the overall optimization procedure. In the 
first step of the procedure, the initial steady-state values 
{Yj, Uj} and the final steady-state values {Yj, Vj], are 
determined, as defined in blocks 226 and 228, respectively. 
In some calculations, both the initial and the final steady- 
state values are required. The initial steady-state values are 
utilized to define the coefficients a 1 , b 1 in a block 228. As 
described above, this utilizes the coefficient scaling of the 
b-coefficients. Similarly, the steady-state values in block 228 
are utilized to define the coefficients a f t b f , it being noted that 
only the b -coefficients are also defined in a block 229. Once 
the beginning and end points are defined, it is then necessary 
to determine the path therebetween. This is provided by 
block 230 for path optimization. There are two methods for 
determining how the dynamic controller traverses this path. 
The first, as described above, is to define the approximate 
dynamic gain over the path from the initial gain to the final 
gain. As noted above, this can incur some instabilities. The 
second method is to define the input values over the horizon 
from the initial value to the final value such that the desired 
value Y 6B f is achieved. Thereafter, the gain can be set for the 
dynamic model by scaling the b-coefificients. As noted 
above, this second method does not necessarily force the 
predicted value of the output y p (t) along a defined path; 
rather, it defines the characteristics of the model as a 
function of the error between the predicted and actual values 
over the horizon from the initial value to the final or desired 
value. This effectively defines the input values for each point 
on the trajectory or, alternatively, the dynamic gain along the 
trajectory. 

[0174] Referring now to FIG. 18a, there is illustrated a 
diagrammatic representation of the manner in which the path 
is mapped through the input and output space. The steady- 
state model is operable to predict both the output steady- 
state value Y^ 1 at a value of k»0, the initial steady-state 
value, and the output steady-state value Y^ 1 at a time t+N 
where k=N, the final steady-state value. At the initial steady- 
state value, there is defined a region 227, which region 227 
comprises a surface in the output space in the proximity of 
the initial steady-state value, which initial steady-state value 
also lies in the output space. This defines the range over 
which the dynamic controller can operate and the range over 
which it is valid. At the final steady-state value, if the gain 
were not changed, the dynamic model would not be valid. 
However, by utilizing the steady-state model to calculate the 
steady-state gain at the final steady-state value and then 
force the gain of the dynamic model to equal that of the 
steady-state model, the dynamic model then becomes valid 



over a region 229, proximate the final steady-state value. 
This is at a value of k«N. The problem that arises is how to 
define the path between the initial and final steady-state 
values. One possibility, as mentioned hereinabove, is to 
utilize the steady-state model to calculate the steady-state 
gain at multiple points along the path between the initial 
steady-state value and the final steady-state value and then 
define the dynamic gain at those points. This could be 
utilized in an optimization routine, which could require a 
large number of calculations. If the computational ability 
were there, this would provide a continuous calculation for 
the dynamic gain along the path traversed between the initial 
steady-state value and the final steady-state value utilizing 
the steady-state gain. However, it is possible that the steady- 
state model is not valid in regions between the initial and 
final steady-state values, i.e., there is a low confidence level 
due to the fact that the training in those regions may not be 
adequate to define the model therein. Therefore, the dynamic 
gain is approximated in these regions, the primary goal 
being to have some adjustment of the dynamic model along 
the path between the initial and the final steady-state values 
during the optimization procedure. This allows the dynamic 
operation of the model to be defined. This is represented by 
a number of surfaces 225 as shown in phantom. 

[0175] Referring now to FIG. 19, there is illustrated a flow 
chart depicting the optimization algorithm. The program is 
initiated at a start block 232 and then proceeds to a function 
block 234 to define the actual input values u a (t) at the 
beginning of the horizon, this typically being the steady- 
state value U B8 . The program then flows to a function block 
235 to generate the predicted values y p (k) over the horizon 
for all k for the fixed input values. The program then flows 
to a function block 236 to generate the error E(k) over the 
horizon for all k for the previously generated y p (k). These 
errors and the predicted values are then accumulated, as 
noted by function block 238. The program then flows to a 
function block 240 to optimize the value of u(t) for each 
value of k in one embodiment. This will result in k-values . 
for u(t). Of course, it is sufficient to utilize less calculations 
than the total k-calculations over the horizon to provide for 
a more efficient algorithm. The results of this optimization 
will provide the predicted change Au(t+k) for each value of 
k in a function block 242. The program then flows to a 
function block 243 wherein the value of u(t+k) for each u 
will be incremented by the value Au(t+k). The program will 
then flow to a decision block 244 to determine if the 
objective function noted above is less than or equal to a 
desired value. If not, the program will flow back along an 
"N" path to the input of function block 235 to again make 
another pass. This operation was described above with 
respect to FIGS. 11a and 116. When the objective function 
is in an acceptable level, the program will flow from 
decision block 244 along the "Y" path to a function block 
245 to set the value of u(t+k) for all u. This defines the path. 
The program then flows to an End block 246. Steady State 
Gain Determination 

[0176] Referring now to FIG. 20, there is illustrated a plot 
of the input space and the error associated therewith. The 
input space is comprised of two variables Xj and x 2 . The 
y-axis represents the function f(x 1( x^. In the plane of xl 
and there is illustrated a region 250, which represents the 
training data set. Areas outside of the region 250 constitute 
regions of no data, i.e., a low confidence level region. The 
function Y will have an error associated therewith. This is 
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represented by a plane 252. However, the error in the plane 
250 is only valid in a region 254, which corresponds to the 
region 250, Areas outside of region 254 on plane 252 have 
an unknown error associated therewith. As a result, when- 
ever the network is operated outside of the region 250 with 
the error region 254, the confidence level in the network is 
low. Of course, the confidence level will not abruptly change 
once outside of the known data regions but, rather, decreases 
as the distance from the known data in the training set 
increases. This is represented in FIG. 21 wherein the con- 
fidence is defined as a(x). It can be seen from FIG. 21 that 
the confidence level a(x) is high in regions overlying the 
region 250. 

[0177] Once the system is operating outside of the training 
data regions, i.e., in a low confidence region, the accuracy of 
the neural net is relatively low. In accordance with one 
aspect of the preferred embodiment, a first principles model 
g(x) is utilized to govern steady-state operation. The switch- 
ing between the neural network model f(x) and the first 
principle models g(x) is not an abrupt switching but, rather, 
it is a mixture of the two. 

[0178] The steady-state gain relationship is defined in 
Equation 7 and is set forth in a more simple manner as 
follows: 



0(3) 



(031) 



[0179] A new output function Y(u) is defined to take into 
account the confidence factor a(u) as follows: 



F(lf)-a(AyfJr))«fJf) 
[0180] where: 

[0181] a(u)=confidence in model f (u) 

[0182] a(u) in the range of 0—1 

[0183] o(u)e{0,l} 

[0184] This will give rise to the relationship: 



(032) 



(033) 



[0185] In calculating the steady-state gain in accordance 
with this Equation utilizing the 110 output relationship Y(u), 
the following will result: 



model. A neural network block 300 is provided for the 
function f(u), a first principle block 302 is provided for the 
function g(u) and a confidence level block 304 for the 
function a(u). The input u(t) is input to each of the blocks 
300-304. The output of block 304 is processed through a 
subtraction block 306 to generate the function l-a(u), 
which is input to a multiplication block 308 for multiplica- 
tion with the output of the first principles block 302. This 
provides the function (l-a(u))*g(u). Additionally, the out- 
put of the confidence block 304 is input to a multiplication 
block 310 for multiplication with the output of the neural 
network block 300. This provides the function f(u)*a(u). The 
output of block 308 and the output of block 310 are input to 
a summation block 312 to provide the output Y(u). 

[0187] Referring now to FIG. 23, there is illustrated an 
alternate embodiment which utilizes discreet switching. The 
output of the first principles block 302 and the neural 
network block 300 are provided and are operable to receive 
the input x(t). The output of the network block 300 and first 
principles block 302 are input to a switch 320, the switch 
320 operable to select either the output of the first principals 
block 302 or the output of the neural network block 300. The 
output of the switch 320 provides the output Y(u). 

[0188] The switch 320 is controlled by a domain analyzer 
322. The domain analyzer 322 is operable to receive the 
input x(t) and determine whether the domain is one that is 
within a valid region of the network 300. If not, the switch 
320 is controlled to utilize the first principles operation in the 
first principles block 302. The domain analyzer 322 utilizes 
the training database 326 to determine the regions in which 
the training data is valid for the network 300. Alternatively, 
the domain analyzer 320 could utilize the confidence factor 
a(u) and compare this with a threshold, below which the first 
principles model 302 would be utilized. 

[0189] Identification of Dynamic Models 

[0190] Gain information, as noted hereinabove, can also 
be utilized in the development of dynamic models. Instead 
of utilizing the user-specified gains, the gains may be 
obtained from a trained steady-stale model. Although 
described hereinabove with reference to Equation 7, a single 
input, single output dynamic model will be defined by a 
similar equation as follows: 



?(0— o^(/-l)-«^(H2>6,«(r-rf-l>4-b 2 u(t-d-2) 



(035) 



[0191] where the dynamic steady-state gain is defined as 
follows: 



Staff)) ,, ^ , _^fl(F(B)) 



(034) 



80) 
d(l-a(u)) 



xg(3) + (l-o(3))- 



6(3) 



[0186] Referring now to FIG. 22, there is illustrated a 
block diagram of the embodiment for realizing the switching 
between the neural network model and the first principles 



1 +Q[ + a-i 



(036) 



[0192] This gain relationship is essentially the same as 
defined hereinabove in Equation 11. Given a time series of 
input and output data, u(t) and y(t), respectively, and the 
steady-state or static gain associated with the average value 
of the input, K'^, the 110 parameters of the dynamic system 
may be defined by minimizing the following cost function: 
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t f (037) 



[0193] where X is a user-specified value. It is noted that the 
second half of Equation 37 constitutes the summation over 
the time series with the value y(t) constituting the actual 
output and the function y p (t) constituting the predicted 
output values. The mean square error of this term is summed 
from an initial time t 4 to a final time t f , constituting the time 
series. The gain value k d basically constitutes the steady- 
state gain of the dynamic model. This optimization is subject 
to the following constraints on dynamic stability: 

0$a 2 <l (038) 
-a 2 -l<a 1 <0 (039) 

[0194] which are conventional constraints. The variable X 
is used to enforce the average steady -state gain, K f 8SJ in the 
identification of the dynamic model. The value K'^ is found 
by calculating the average value of the steady-state gain 
associated with the neural network over the time horizon tj 
to t f . Given the input time series u(tj) to u(tf), the }C 6S is 
defined as follows: 




(040) 



[0195] For a large value of X, the gain of the steady-state 
and dynamic models are forced to be equal. For a small 
value of X, the gain of the dynamic model is found inde- 
pendently from that of the steady-state model. For X=0, the 
optimization problem is reduced to a technique commonly 
utilized in identification of output equation-based models, as 
defined in L. Ljung, "System Identification: Theory for the 
User," Prentice-Hall, Englewood Clifis, N.J. 1987. 

[0196] In defining the dynamic model in accordance with 
Equation No. 37, it is recognized that only three parameters 
need to be optimized, the a a parameter, the a 2 parameter and 
the ratio of b x and b 2 . This is to be compared with the 
embodiment described hereinabove with reference to FIG. 
2, wherein the dynamic gain was forced to be equal to the 
steady-state gain of the static model 20. By utilizing the 
weighting factor X and minimizing the cost function in 
accordance with Equation 37 without requiring the dynamic 
gain k d to equal the steady-state gain K' ss of the neural 
network, some latitude is provided in identifying the 
dynamic model. 

[0197] In the embodiment described above with respect to 
FIG. 2, the model was identified by varying the b-values 
with the dynamic gain forced to be equal to the steady-state 
gain. In the embodiment illustrated above with respect to 
Equation 37, the dynamic gain does not necessarily have to 
equal the steady-state gain K'„, depending upon the value of 
X defining the weighting factor. 

[0198] The above noted technique of Equation 37 provides 
for determining the a's and b's of the dynamic model as a 
method of identification in a particular localized region of 



the input space. Once the a's and b's of the dynamic model 
are known, this determines the dynamics of the system with 
the only variation over the input space from the localized 
region in which the dynamic step test data was taken being 
the dynamic gain k d . If this gain is set to a value of one, then 
the only component remaining are the dynamics. Therefore, 
the dynamic model, once defined, then has its gain scaled to 
a value of one, which merely requires adjusting the b-values. 
This will be described hereinbelow. After identification of 
the model, it is utilized as noted hereinabove with respect to 
the embodiment of FIG. 2 and the dynamic gain can then be 
defined for each region utilizing the static gain. 

[0199] Steady-State Model Identification 

[0200] As noted hereinabove, to optimize and control any 
process, a model of that process is needed. The present 
system relies on a combination of steady-state and dynamic 
models. The quality of the model determines the overall 
quality of the final control of the plant. Various techniques 
for training a steady-state model will be described herein- 
below. 

[0201] Prior to discussing the specific model identification 
method utilized in the present embodiment, it is necessary to 
define a quasi-steady-state model. For comparison, a steady- 
state model will be defined as follows: 

[0202] Steady-State Models: 

[0203] A steady-state model is represented by the static 
mapping from the input, u(t) to the output y(t) as defined by 
the following equation: 

?(0-W)) (041) 

[0204] where F(u(t)) models the steady-state map- 
ping of the process and u(t) 6 R m and y(t) E R n 
represent the inputs and outputs of a given process. 
It should be noted that the input, u(t), and the output, 
y(t), are not a function of time and, therefore, the 
steady-state mapping is independent of time. The 
gain of the process must be defined with respect to a 
point in the input space. Given the point u(t), the gain 
of process is defined as: 



dy I (042) 
cm) = -£: 

du Uj 

[0205] where G is a R mxn matrix. This gain is 

equivalent to the sensitivity of the function F(u(t)). 

[0206] Quasi-Stcady-Statc Models: 

[0207] A steady-state model by definition contains no time 
information. In some cases, to identify steady-state models, 
it is necessary to introduce time information into the static 
model: 

f(fy-G<#M) (043) 

[0208] where: 

*(f,dH"iit-dMf-dJ ■ • • "cC-O] (044) 
[0209] The variable d { represents the delay associated 
with the i th input. In the quasi-steady-state model, the 
static mapping of G(u(t)) is essentially equal to the 
steady-state mapping F(u(t)). The response of such a 
model is illustrated in FIG. 24. In FIG. 24, there is 
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illustrated a single input Uj(t) and a single output 
yj(t), this being a single output, single input system. 
There is a delay, or dead time d, noted between the 
input and the output which represents a quasi -steady- 
state dynamic. It is noted, however, that each point 
on the input u 1 (t) corresponds to a given point on the 
output yj(t) by some delay d. However, when ul(t) 
makes a change from an initial value to a final value, 
the output makes basically an instantaneous change 
and follows it. With respect to the quasi-steady-state 
model, the only dynamics that are present in this 
model is the delay component. 

[0210] Identification of Delays in Quasi -Steady-State 
Models 

[0211] Given data generated by a quasi -steady-state model 
of the form 

y{t)=Gm-d)), (045) 

[0212] where d is the dead-time or delay noted in FIG. 24, 
the generating function G( ) is approximated via a neural 
network training algorithm (nonlinear regression) when d is 
known. That is, a function G( ) is fitted to a set of data points 
generated from G( ), where each data point is a u(t), y(t) pair. 
The present system concerns time-scries data, and thus the 
dataset is indexed by t. The data set is denoted by D. 

[0213] In process modeling, exact values for d are ordi- 
narily not critical to the quality of the model; approximate 
values typically suffice. Prior art systems specified a method 
for approximating d by training a model with multiple 
delays per input, and picking the delay which has the largest 
sensitivity (average absolute partial derivative of output 
w.r.t. input). In these prior art systems, the sensitivity was 
typically determined by manipulating a given input and 
determining the effect thereof on the output. By varying the 
delay, i.e., taking a different point of data in time with 
respect to a given y(t) value, a measure of sensitivity of the 
output on the input can be determined as a function of the 
delay. By taking the delay which exhibits the largest sensi- 
tivity, the delay of the system can be determined. 

[0214] The disadvantage to the sensitivity technique is that 
it requires a number of passes through the network during 
training in order to determine the delay. This is an iterative 
technique. In accordance with the present system, the 
method for approximating the delay is done utilizing a 
statistical method for examining the data, as will be 
described in more detail hereinbelow. This method is per- 
formed without requiring actual neural network training 
during the determination operation. The method of the 
present embodiment examines each input variable against a 
given output variable, independently of the other input 
variables. Given d ; for an input variable u., the method 
measures the strength of the relationship between u^t-di) 
and y(t). The method is fast, such that many di values may 
be tested for each u r 

[0215] The user supplies d i mIn and di )inflX values for each 
u t . The strength of the relationship between u^t-dj and y(t) 
is computed for each d { between d i min and d i?nQX (inclusive). 
The dj yielding the strongest relationship between u^t-d^) 
and y(t) is chosen as the approximation of the dead-time for 
that input variable on the given output variable. The strength 
of the relationship between u^t-dj and y(t) is defined as the 
degree of statistical dependence between u^t-d;) and y(t). 



The degree of statistical dependence between u^t-dj and 
y(t) is the degree to which u^t-dj and y(t) are not statisti- 
cally independent. 

[0216] Statistical dependence is a general concept. As long 
as there is any relationship whatsoever between two vari- 
ables, of whatever form, linear or nonlinear, the definition of 
statistical independence for those two variables will fail. 
Statistical independence between two variables x a (t) and 
x 2 (0 is defined as: 



(046) 



[0217] where p(x 2 (t)) is the marginal probability density 
function of Xj(t) and p(x 1 (t),x 2 (/)) is the joint probability 
density function (x^O-u^t-aV) and x 2 =y(t)); that is, the 
product of the marginal probabilities is equal to the joint 
probability. If they are equal, this constitutes statistical 
independence, and the level of inequality provides a mea- 
sure of statistical dependence. 

[0218] Any measure £(x 1 (t)>x 2 (t)) which has the following 
property ("Property 1") is a suitable measure of statistical 
dependence: 

[0219] Property 1: f(Xi(t),x 2 (/)) is 0 if and only if 
Equation 46 holds at each data point, and f>0 oth- 
erwise. In addition, the magnitude of f measures the 
degree of violation of Equation 46 summed over all 
data points. 

[0220] Mutual information (MI) is one such measure, and 
is defined as: 



mi = £p<*i(0,*2(0)log( 



P(*l(0,*2(0) \ 

p{x x {t))p{x 2 {t))) 



(047) 



[0221] Property 1 holds for MI. Theoretically, there is no 
fixed maximum value of MI, as it depends upon the distri- 
butions of the variables in question. As explained hereinbe- 
low, the maximum, as a practical matter, also depends upon 
the method employed to estimate probabilities. Regardless, 
MI values for different pairs of variables may be ranked by 
relative magnitude to determine the relative strength of the 
relationship between each pair of variables. Any other 
measure f having Property 1 would be equally applicable, 
such as the sum of the squares of the product of the two sides 
of Equation 58: 



ssd = KpUi (0. *2(0) - P(xi(0)p(x 2 (t)))\ 2 



(048) 



[0222] However, MI (Equation 47) is the preferred method 
in the disclosed embodiment. 

[0223] Statistical Dependence vs. Correlation 

[0224] For purposes of the present embodiment, the 
method described above, i.e., measuring statistical depen- 
dence, is superior to using linear correlation. The definition 
of linear correlation is well-known and is not slated herein. 
Correlation ranges in value from -1 to 1, and its magnitude 
indicates the degree of linear relationship between variables 
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Xj(t) and x 2 (t). Nonlinear relationships are not detected by 
correlation. For example, y(t) is totally determined by x(t) in 
the relation 

y(0'A0- (049) 

[0225] Yet, if x(t) varies symmetrically about zero, then: 

corr(y(t)rtt))-0. (050) 

[0226] That is, correlation detects no linear relationship 
because the relationship is entirely nonlinear. Conversely, 
statistical dependence registers a relationship of any kind, 
linear or nonlinear, between variables. In this example, 
MI(y(t),x(t)) would calculate to be a large number. 

[0227] Estimation Probabilities 

[0228] An issue in computing MI is how to estimate the 
probability distributions given a dataset D. Possible methods 
include kernel estimation methods, and binning methods. 
The preferred method is a binning method, as binning 
methods are significantly cheaper to compute than kernel 
estimation methods. 

[0229] Of the binning techniques, a very popular method 
is that disclosed in A. M. Fraser and Harry L. Swinney, 
"Independent Coordinates for Strange Attractors in Mutual 
Information,"P/jysica/ Review A, 33(2): 1134-1140, 1986. 
This method makes use of a recursive quadrant-division 
process. 

[0230] The present method uses a binning method whose 
performance is highly superior to that of the Fraser method. 
The binning method used herein simply divides each of the 
two dimensions (u £ (t-dj) and y(t)) into a fixed number of 
divisions, where N is a parameter which may be supplied by 
the user, or which defaults to sqrt(#datapoints/20). The 
width of each division is variable, such that an (approxi- 
mately) equal number of points fall into each division of the 
dimension. Thus, the process of dividing each dimension is 
independent of the other dimension. 

[0231] In order to implement the binning procedure, it is 
first necessary to define a grid of data points for each input 
value at a given delay. Each input value will be represented 
by a time series and will therefore be a series of values. For 
example, if the input value were u 1 (t), there would be a time 
series of these u a (t) values, u^tj, u^tj . . . u^tf). There 
would be a time series u a (t) for each output value y(t). For 
the purposes of the illustration herein, there will be consid- 
ered only a single output from y(t), although it should be 
understood that a multiple input, multiple output system 
could be utilized. 

[0232] Referring now to FIG. 25, there is illustrated a 
diagrammatic view of a binning method. In this method, a 
single point generated for each value of U;(t) for the single 
value y(t). All of the data in the time series Uj(t) is plotted in 
a single grid. This time series is then delayed by the delay 
value d. to provide a delay value u^t-d.). For each value of 
j from d j min to dj jmax , there will be a grid generated. There 
will then be a mutual information value generated for each 
grid to show the strength of the relationship between that 
particular delay value dj and the output y(t). By continually 
changing the delay dj for the time series u { (t), a different MI 
value can be generated. 

[0233] In the illustration of FIG. 25, there are illustrated 
a plurality of rows and a plurality columns with the data 



points disposed therein with the x-axis labeled u a (t-dj) and 
the y-axis labeled y(t). For a given column 352 and a given 
row 354, there is defined a single bin 350. As described 
above, the grid lines are variable such that the number of 
points in any one division is variable, as described herein- 
above. Once the grid is populated, then it is necessary to 
determine the MI value. This MI value for the binning grid 
or a given value of dj is defined as follows: 



^"Z£^'^o))iog ^ (0)pCn(/)) 



(051) 



[0234] where p(x 1 (I),x 2 (j)) is equal to the number of 
points in a particular bin over the total number of points in 
the grid, p(x a (I)) & equal to the number of points in a column 
over the total number of points and p(x 2 (j)) is equal to the 
number of data points in a row over the total number of data 
points in the grid and n is equal to the number of rows and 
M is equal to the number of columns. Therefore, it can be 
seen that if the data was equally distributed around the grid, 
the value of MI would be equal to zero. As the strength of 
the relationship increases as a function of the delay value, 
then it would be noted that the points tend to come together 
in a strong relationship, and the value of MI increases. The 
delay dj having the strongest relationship will therefore be 
selected as the proper delay for that given u^t), 

[0235] Referring now to FIG. 26, there is illustrated a 
block diagram depicting the use of the statistical analysis 
approach. The statistical analysis is defined in a block 353 
which receive both the values of y(t) and u(t). This statistical 
analysis is utilized to select for each u^t) the appropriate 
delay dj. This, of course, is for each y(t). The output of this 
is stored in a delay register 355. During training of a 
non-linear neural network 357, a delay block 359 is provided 
for selecting from the data set of u(t) for given u^t) an 
appropriate delay and introducing that delay into the value 
before inputting it to a training block 358 for training the 
neural network 357. The training block 358 also utilizes the 
data set for y(t) as target data. Again, the particular delay for 
the purpose of training is defined by statistical analysis block 
353, in accordance with the algorithms described herein- 
above. 

[0236] Referring now to FIG. 27, there is illustrated a flow 
chart depicting the binning operation. The procedure is 
initiated at a block 356 and proceeds to a block 360 to select 
a given one of the outputs y(t) for a multi-output system and 
then to a block 362 to select one of the input values u £ (l— dj). 
It then flows to a function block 363 to set the value of dj to 
the minimum value and then to a block 364 to perform the 
binning operation wherein all the points for that particular 
delay dj are placed onto the grid. The MI value is then 
calculated for this grid, as indicated by a block 365. The 
program then proceeds to a decision block 366 to determine 
if the value of dj is equal to the maximum value dj max . If not, 
this value is incremented by a block 367 and then proceeds 
back to the input of block 364 to increment the next delay 
value for u^t-dj). This continues until the delay has varied 
from d i min through dj max . The program then flows to the 
decision block 368 to determine if there are additional input 
variables. If so, the program flows to a block 369 to select 
the next variable and then back to the input of block 363. If 
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not, the program flows to a block 370 to select the next value 
of y(t). This will then flow back to the input of function 
block 360 until all input variables and output variables have 
been processed. The program will then flow to an END 
block 371. 

[0237] Identification of Steady-State Models Using Gain 
Constraints: 

[0238] In most processes, bounds upon the steady-state 
gain are known either from the first principles or from 
practical experience. Once it is assumed that the gain 
information is known, a method for utilizing this knowledge 
of empirically-based models will be described herein. If one 
considers a parameterized quasi-steady-state model of the 
form: 

?(/)-SC*,irM)) (052) 
[0239] where w is a vector of free parameters (typically 
referred to as the weights of a neural network) and N(w, 
u(t-d)) represents a continuous function of both w and 
u(t-d). A feedforward neural network as described herein- 
above represents an example of the nonlinear function. A 
common technique for identifying the free parameters w is 
to establish some type of cost function and then minimize 
this cost function using a variety of different optimization 
techniques, including such techniques as steepest descent or 
conjugate gradients. As an example, during training of 
feedforward neural networks utilizing a backpropogation 
algorithm, it is common to minimize the mean squared error 
over a training set, 

p (053) 



[0240] where P is the number of training patterns, y d (t) is 
the training data or target data, y(t) is the predicted output 
and J(w) is the error. 

[0241] Constraints upon the gains of steady-state models 
may be taken into account in determining w by modifying 
the optimization problem. As noted above, w is determined 
by establishing a cost function and then utilizing an opti- 
mization technique to minimize the cost function. Gain 
constraints may be introduced into the problem by specify- 
ing them as part of the optimization problem. Thus, the 
optimization problem may be reformulated as: 



min(J(tf)) (054) 

[0242] subject to 

c l m)xcm))<G h m)) (° 55 > 

C,fj3(2)) < C(b(2)) < G h &{2)) (056) 

i (057) 

CjOT < G{u(P)) < G h (u{P}) (058) 



[0243] where G^t)) is the matrix of the user-specified 
lower gain constraints and G h (u(t)) are the upper gain 
constraints. Each of the gain constraints represents the 
enforcement of a lower and upper gain on a single one of the 



input -output pairs of the training set, i.e., the gain is bounded 
for each input-output pair and can have a different value. 
These are what are referred to as "hard constraints." This 
optimization problem may be solved utilizing a non-linear 
programming technique. 

[0244] Another approach to adding the constraints to the 
optimization problem is to modify the cost function, i.e., 
utilize some type of soft constraints. For example, the 
squared error cost function of Equation 53 may be modified 
to account for the gain constraints in the gain as follows: 



A 2 (059) 

Aw) = 2j (?(0 -?,(/)) + 

p 

A£ (tf (Gi(3t/)) - G(u(t))) + H(G(H(t)) - CmOHO))) 



[0245] where H( # ) represents a non-negative penalty func- 
tion for violating the constraints and X is a user-specified 
parameter for weighting the penalty. For large values of X, 
the resulting model will observe the constraints upon the 
gain. In addition, extra data points which are utilized only in 
the second part of the cost function may be added to the 
historical data set to effectively fill voids in the input space. 
By adding these additional points, proper gain extrapolation 
or interpolation can be guaranteed. In the preferred embodi- 
ment, the gain constraints are held constant over the entire 
input space. 

[0246] By modifying the optimization problem with the 
gain constraints, models that observe gain constraints can be 
effectively trained. By guaranteeing the proper gain, users 
will have greater confidence that an optimization and control 
system based upon such a model will work properly under 
all conditions. 

[0247] One prior art example of guaranteeing global posi- 
tive or negative gain (monotonicity) in a neural network is 
described in J. Sill & Y S, Abu-Mostafa, "Monotonicity 
Hints," Neural Information Processing Systems, 1996. The 
technique disclosed in this reference relies on adding an 
additional term to the cost function. However, this approach 
can only be utilized to bound the gain to be globally positive 
or negative and is not used to globally bound the range of the 
gain, nor is it utilized to locally bound the gain (depending 
on the value of u(t)). 

[0248] Identification of SS Model with Dynamic Data 

[0249] Referring now to FIG. 28, there is illustrated a 
block diagram of a prior art Weiner model, described in M. 
A. Henson and D. F. Seborg, "Nonlinear Process Control," 
Prentice Hall PTR, 1997, Chapter 2, ppll-110. In the Weiner 
model, a non-linear model 376 is generated. This non-linear 
model is a steady-state model. This steady-state model may 
be trained on input data u(t) to provide the function y(t)« 
f(u(t)) such that this is a general non-linear model. However, 
the input u(t) is processed through a linear dynamic model 
374 of the system, which linear dynamic model 374 has 
associated therewith the dynamics of the system. This pro- 
vides on the output thereof a filtered output u(t) which has 
the dynamics of the system impressed thereupon. This 
constitutes the input to the non-linear model 376 to provide 
on the output a prediction y(t). 
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[0250] Referring now to FIG. 29, there is illustrated a 
block diagram of the training method of the present embodi- 
ment. A plant 378 is provided which can represent any type 
of system to be modeled, such as a boiler or a chemical 
process. There are various inputs provided to the plant in the 
form of u(t). This will provide an actual output y^t). 
Although not illustrated, the plant has a number of measur- 
able state variables which constitute the output of various 
sensors such as flow meters, temperature sensors, etc. These 
can provide data that is utilized for various training opera- 
tions, these state outputs not illustrated, it being understood 
that the outputs from these devices can be a part of the input 
training data set. 

[0251] A steady-state neural network 379 is provided 
which is a non-linear network that is trained to represent the 
plant. A neural network typically contains an input layer and 
an output layer and one or more hidden layers. The hidden 
layers provide the mapping for the input layers to the output 
layers and provide storage for the stored representation of 
the plant. As noted hereinabove, with a sufficient amount of 
steady-state data, an accurate steady-state model can be 
obtained. However, in a situation wherein there is very little 
steady-state data available, the accuracy of a steady-state 
model with conventional training techniques is questionable. 
As will be described in more detail hereinbelow, the training 
method of the present embodiment allows training of the 
neural network 374, or any other empirical modeling 
method to learn the steady-state process model from data 
that has no steady-state behavior, i.e., there is a significant 
dynamic component to all training data. 

[0252] Typically, a plant during operation thereof will 
generate historical data. This historical data is collected and 
utilized to later train a network. If there is little steady-state 
behavior exhibited in the input data, the present embodiment 
allows for training of the steady -state model. The input data 
u(t) is input to a falter 381 which is operable to impress upon 
the input data the dynamics of the plant 378 and the training 
data set. This provides a filtered output u f (t) which is input 
to a switch 380 for input to the plant 378. The switch 380 is 
operable to input the unfiltered input data u(t) during opera- 
tion of the plant, or the filtered input data u f (t) during 
training into the neural network 379. As will be described 
hereinbelow, the u(t) input data, prior to being filtered, is 
generated as a separate set of dynamic training data by a step 
process which comprises collecting step data in a local 
region. The filter 381 has associated therewith a set of 
system dynamics in a block 382 which allows the filter 381 
to impress the dynamics of the system onto the input training 
data set. Therefore, during training of the neural network 
379, the filtered data u(t) is utilized to train the network such 
that the neural network 379 provides an output y(t) which is 
a function of the filtered data -tf(t) or: 

y p (')=/ttt f (0) (060) 
[0253] Referring now to FIG. 30, there is illustrated a 
diagrammatic view of the training data and the output data. 
The training data is the actual set of training data which 
comprises the historical data. This is the u(t) data which can 
be seen to vary from point to point. The problem with some 
input data in a training set of data, if not all data, is that it 
changes from one point to another and, before the system 
has "settled," it will change again. That is, the average time 
between movements in u(t) is smaller than T„, where T„ is 
the time for y(t) to reach steady-state. As such, the corre- 



sponding output data y(t) will constitute dynamic data or 
will have a large dynamic component associated therewith. 
In general, the presence of this dynamic information in the 
output data must be accounted for to successfully remove 
the dynamic component of the data and retain the steady- 
state component of the steady-state neural network 379. 

[0254] As will be described in more detail hereinbelow, 
the present embodiment utilizes a technique whereby the 
actual dynamics of the system which are inherent in the 
output data y(t) are impressed upon the input data u(t) to 
provide filtered input data u(t). This data is scaled to have a 
gain of one, and the steady-state model is then trained upon 
this filtered data. As will also be described in more detail 
hereinbelow, the use of this filtered data essentially removes 
the dynamic component from the data with only the steady- 
state component remaining. Therefore, a steady-state model 
can be generated. 

[0255] Referring now to FIG. 31, there is illustrated a flow 
chart depicting the training procedure for the neural network 
379 of FIG. 29 for a single output. As noted above, the 
neural network is a conventional neural network comprised 
of an input layer for receiving a plurality of input vectors, an 
output layer for providing select predicted outputs, and one 
or more hidden layers which are operable to map the input 
layer to the output layer through a stored representation of 
the plant 378. This is a non-linear network and it is trained 
utilizing a training data set of input values and target output 
values. This is, as described hereinabove, an iterative pro- 
cedure utilizing algorithms such as the backpropagation 
training technique. Typically, an input value is input to the 
network during the training procedure, and also target data 
is provided on the output. The results of processing the input 
data through the network are compared to the target data, 
and then an error generated. This error, with the backpropa- 
gation technique, is then back propagated through the net- 
work from the output to the input to adjust the weights 
therein, and then the input data then again processed through 
the network and the output compared with the target data to 
generate a new error and then the algorithm readjusts the 
weights until they are reduced to an acceptable level. This 
can then be used for all of the training data with multiple 
passes required to minimize the error to an acceptable level, 
resulting in a trained network that provides therein a stored 
representation of the plant. 

[0256] As noted above, one of the disadvantages to con- 
ventional training methods is that the network 379 is trained 
on the set of historical input data that can be incomplete, or 
have some error associated therewith. The incompleteness 
of the historical data may result in areas in the input space 
on which the network is not trained. The network, however, 
will extrapolate its training data set during the training 
operation, and actually provide a stored representation 
within that portion of the input space in which data did not 
exist. As such, whenever input data is input to the network 
in an area of the input space in which historical input data 
did not exist during training, the network will provide a 
predicted output value. This, however, effectively decreases 
the confidence level in the result in this region. Of course, 
whenever input data is input to the network in a region that 
was heavily populated with input data, the confidence level 
is relatively high. 

[0257] Another source of error, as noted hereinabove, is 
the dynamic component of the data. If the historical data that 
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forms the training data set is dynamic in nature, i.e., it is 
changing in such a manner that the output never settles down 
to a steady-state value, this can create some errors when 
utilizing this data for training. The reason for this is that the 
fundamental assumption in training a steady-state neural 
network with an input training data set is that the data is 
steady-state data. The training procedure of the present 
embodiment removes this error. 

[0258] Referring further to FIG. 31, the flow chart is 
initiated at a block 384 and then proceeds to a block 386 
wherein dynamic data for the system is collected. In the 
preferred embodiment, this is in the form of step test data 
wherein the input is stepped between an initial value and a 
final value multiple times and output data taken from the 
plant under these conditions. This output data will be rich in 
dynamic content for a local region of the input space. An 
alternative method is to examine the historical data taken 
during the operation of the plant and examine the data for 
movements in the manipulated variables (MVs) or the 
dynamic variables (DVs). These variables are then utilized 
for the purpose of identifying the dynamic model. However, 
the preferred model is to utilize a known input that will 
result in the dynamic change in the output. Of course, if 
there are no dynamics present in the output, then this will 
merely appear as a steady-state value, and the dynamic 
model will have filter values of a«0 and b-0. This will be 
described hereinbelow. 

[0259] The step test data, as will be described hereinbe- 
low, is data that is taken about a relatively small region of 
the input space. This is due to the fact that the variables are 
only manipulated between two values, and initial steady- 
state value and a final value, in a certain region of the input 
space, and the data is not taken over many areas of the input 
space. Therefore, any training set generated will represent 
only a small portion of the input space. This will be 
described in more detail hereinbelow. It should be noted that 
these dynamics in this relatively small region of the input 
space will be utilized to represent the dynamics over the 
entire input space. A fundamental presumption is that the 
dynamics at any given region remain substantially constant 
over the entire input space with the exception of the dynamic 
gain varying. 

[0260] Once the dynamic data has been collected for the 
purpose of training, this dynamic training data set is utilized 
to identify the dynamic model of the system. If, of course, 
a complete steady-state data set were available, there would 
be a reduced need for the present embodiment, although it 
could be utilized for the purpose of identifying the dynamics 
of the system. The flow chart then proceeds to a block 387 
wherein the dynamics of the plant are identified. In essence, 
a conventional model identification technique is utilized 
which models the dynamics of the plant. This is a linear 
model which is defined by the following equation: 

[0261] In the above-noted model of Equation 61, the 
values of a 1( a 2 , b 2 and b 2 define the parameters of the model 
and are defined by training this model. This operation will be 
described in detail hereinbelow; however, once trained, this 
model will define the dynamic model of the plant 378 as 
defined by the dynamics associated with the dynamic train- 
ing data set at the location in the input space at which the 
data was taken. This will, of course, have associated there- 



with a dynamic gain, which dynamic gain will change at 
different areas in the input space. 

[0262] Once the dynamic model has been identified uti- 
lizing the dynamic training data set, i.e., the a's and b's of 
the model have been determined, the program will flow to a 
function block 388 to determine the properties of a dynamic 
pre-filter model, which is operable to process the input 
values u(t) through the dynamic model to provide a filtered 
output u (t) on the output which is, in effect, referred to as 
a "filtered" input in accordance with the following equation: 

Tt f 0>-ailj('-l)-a 2 uC'-l)+&i«(0+M('-3) (062) 

[0263] wherein the values of a a and a 2 are the same as in 
the dynamic model of the plant, and the values of b 1 and 5 2 
are adjusted to set the gain to a value of zero. 

[0264] The pre-filter operation is scaled such that the gain 
of the dynamic model utilized for the pre-filter operation is 
set equal to unity. The b-values are adjusted to provide this 
gain scaling operation. The gain is scaled in accordance with 
the following: 



na = 1 = 

1 + Q\ +fl2 



[0265] If the gain were not scaled, this would require some 
adjustment to the steady-state model after training of the 
steady-state model. For example, if the gain of the model 
were equal to "two," this would require that the steady-state 
model have a gain adjustment of "one-half after training. 

[0266] After the filter values have been determined, i.e., 
the b-values with the gain set equal to one, then the input 
values u(t) for the historical data are processed through the 
pre-filter with the gain set equal to one to yield the value of 
u(t), as indicated by a function block 390. At this point, the 
dynamics of the system are now impressed upon the his- 
torical input data set, i.e., the steady-state component has 
been removed from the input values. These input values u(t) 
are now input to the neural network in a training operation 
wherein the neural network 378 is trained upon the filtered 
input values over the entire input space (or whatever portion 
is covered by the historical data). This data u(t) has the 
dynamics of the system impressed thereupon, as indicated 
by block 391. The significance of this is that the dynamics 
of the system have now been impressed upon the historical 
input data and thus removed from the output such that the 
only thing remaining is the steady-state component. There- 
fore, when the neural network 378 is trained on the filtered 
output, the steady-state values are all that remain, and a valid 
steady-state model is achieved for the neural network 378. 
This steady-state neural network is achieved utilizing data 
that has very little steady-state nature. Once trained, the 
weights of the neural network are then fixed, as indicated by 
a function block 392, and then the program proceeds to an 
END block 394. 

[0267] Referring now to FIG. 32, there is illustrated a 
diagrammatic view of the step test. The input values u(t) are 
subjected to a step response such that they go from an initial 
steady-state value u^t) to a final value Uf(t). This results in 
a response on the output y(t), which is rich in dynamic 
content. The step test is performed such that the value of u(t) 
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is increased from u £ (t) to u^t), and then decreased back to 
Uj(t), preferably before the steady -state value has been 
reached. The dynamic model can then be identified utilizing 
the values of u(t) presented to the system in the step test and 
the output values of y(t), these representing the dynamic 
training data set. This information is utilized to identify the 
model, and then the model utilized to obtain the pre-filtered 
values of u(t) by passing the historical input data u(t) 
through the identified model with the dynamic gain of the 
model set equal to one. Again, as noted above, the values of 
b are adjusted in the model such that the gain is set equal to 
one. 

[0268] Referring now to FIG. 33, there is illustrated a 
diagrammatic view of the relationship between the u(t) and 
the u(t), indicating that the gain is set equal to one. By setting 
the gain equal to one, then only the dynamics determined al 
the "training region" will be impressed upon the historical 
input data which exists over the entire input space. If the 
assumption is true that the only difference between the 
dynamics between given regions and the input space is the 
dynamic gain, then by setting the gain equal to one, the 
dynamics at the given region will be true for every other 
region in the input space. 

[0269] Referring now to FIG. 34, there is illustrated a 
block diagram of the system for training a given output. The 
training system for the embodiment described herein with 
respect to impressing the dynamics of the plant onto the 
historical input data basically operates on a single output. 
Most networks have multiple inputs and multiple outputs 
and are referred to as MI MO (multi-input, multi -output) 
networks. However, each output will have specific dynamics 
that are a function of the inputs. Therefore, each output must 
have a specific dynamic model which defines the dynamics 
of that output as a function of the input data. Therefore, for 
a given neural network 400, a unique pre-filter 402, will be 
required which receives the input data on the plurality of 
input lines 404. This pre-filter is operable to incorporate a 
dynamic model of the output y(t) on the input u(t). This will 
be defined as the function: 

[0270] This represents the dynamic relationship between 
the inputs and a single output with the gain set equal to unity. 

y(0-/ d ,oc,U (u(0) (064) 

[0271] Referring now to FIG. 35, there is illustrated a 
dynamic representation of a MIMO network being modeled. 
In a MIMO network, there will be required a plurality of 
steady-state neural networks 410, labeled NN lt NN 2 , . . . 
NN M . Each one is associated with a separate output y^t), 
v 2(0> • • ■ YmO)- F° r eacn °f me neural networks 410, there 
will be a pre-filter or dynamic model 412 labeled Dyn^ 
Dyn 2 . . . Dyn M . Each of these models 412 receives on the 
input thereof the input values u(t), which constitutes all of 
the inputs u^t), u 2 (t), . . . u^t). For each of the neural 
networks 410, during the training operation, there will also 
be provided the dynamic relationship between the output 
and the input u(t) in a block 14. This dynamic relationship 
represents only the dynamic relationship between the asso- 
ciated one of the outputs y a (t), y 2 (r), . . . y M (0* Therefore, 
each of these neural networks 410 can be trained for the 
given output. 

[0272] Referring now to FIG. 36, there is illustrated the 
block diagram of the predicted network after training. In this 
mode, all of the neural networks 410 will now be trained 



utilizing the above-noted method of FIG. 30, and they will 
be combined such that the input vector u(t) will be input to 
each of the neural networks 410 with the output of each of 
the neural networks comprising one of the outputs y 1 (t), 

y 2 (0. • • • y M «. 

[0273] Graphical Interface for Model Identification 

[0274] Referring now to FIG. 37, there is illustrated a 
graphical user interface (GUI) for allowing the user to 
manipulate data on the screen, which manipulated data then 
defines the parameters of the model identification procedure. 
In FIG. 39 there are illustrated two input values, one labeled 
"reflux" and one labeled "steam." The reflux data is repre- 
sented by a curve 400 whereas the steam data is represented 
by curve 402. These curves constitute dynamic step data 
which would have corresponding responses on the various 
outputs (not shown). As noted hereinabove with respect to 
FIG. 34, the output would have a dynamic response as a 
result of the step response of the input. 

[0275] The user is presented the input data taken as a result 
of the step test on the plant and then allowed to identify the 
model from this data. The user is provided a mouse or 
similar pointing device (not shown) to allow a portion of one 
or more of the data values to be selected over a user defined 
range. In FIG. 37, there is illustrated a box 404 in phantom 
about a portion of the input reflux data which is generated by 
the user with the pointing device. There is also illustrated a 
box 406 in phantom about a portion of the steam data on 
curve 402. The portion of each of the curves 400 and 402 
that is enclosed within the respective boxes 404 and 406 is 
illustrated in thick lines as "selected" data. As noted here- 
inabove, the step test data is taken in a particular localized 
portion of the input space, wherein the input space is defined 
by the various input values in the range over which the data 
extends. By allowing the user the versatility of selecting 
which input data is to be utilized for the purpose of identi- 
fying the model, the user is now permitted the ability to 
manipulate the input space. Once the user has selected the 
data that is to be utilized and the range of data, the user then 
selects a graphical button 410 which will then perform an 
"identify" operation of the dynamic model utilizing the 
selected information. 

[0276] Referring now to FIG. 38, there is illustrated a 
flowchart depicting the general identify operation described 
above with respect to FIG. 39. The program is initiated at a 
function block 412 which indicates a given data set, which 
data set contains the step test data for each input and each 
output. The program then flows to a function block 414 to 
display the step test data and then to a function block 416 to 
depict the operation of FIG. 38 wherein the user graphically 
selects portions of the display step test data. The program 
then flows to a function block 418 wherein the data set is 
modified with the selected step test data. It is necessary to 
modify the data set prior to performing the identification 
operation, as the identification operation utilizes the avail- 
able data set. Therefore, the original data set will be modi- 
fied to only utilize the selected data, i.e., to provide a 
modified data set for identification purposes. The program 
will then flow to a function block 420 wherein the model 
will be identified utilizing the modified data set. 

[0277] Referring now to FIG. 39, there is illustrated a 
second type of graphical interface. After a model has been 
created and identified, it is then desirable to implement the 
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model in a control environment to control the plant by 
generating new input variables to change the operation of 
the plant for the purpose of providing a new and desired 
output. However, prior to placing the model in a "run-time" 
mode, it may be desirable to run a simulation on the model 
prior to actually incorporating it into the run time mode. 
Additionally, it may be desirable to graphically view the 
system when running to determine how the plant is operat- 
ing in view of what the predicted operation is and how that 
operation will go forward in the future. 

[0278] Referring further to FIG. 40, there is illustrated a 
plot of a single manipulatable variable (MV) labeled 
"reflux." The plot illustrates two sections, a historical sec- 
tion 450 and a predictive section 452. There is illustrated a 
current time line 454 which represents the actual value at 
that time. The x-axis is comprised of the steps and horizontal 
axis illustrates the values for the MV. In this example of 
FIG. 41, the system was initially disposed at a value of 
approximately 70.0. In historical section 450, a bold line 456 
illustrates the actual values of the system. These actual 
values can be obtained in two ways. In a simulation mode, 
use the actual values that are input to the simulation model. 
In the run-time mode, they are the set point values input to 
the plant. 

[0279] In a dynamic system, any change in the input will 
be made in a certain manner, as described hereinabove. It 
could be a step response or it could have a defined trajectory. 
In the example of FIG. 40, the change in the MV is defined 
along a desired trajectory 458 wherein the actual values are 
defined along a trajectory 460. The trajectory 458 is defined 
as the "set point" values. In a plant, the actual trajectory may 
not track the setpoints (desired MVs) due to physical limi- 
tations of the input device, time constraints, etc. In the 
simulation mode, these are the same curve. It can be seen 
that the trajectory 460 continues up to the current time line 
454. After the current time line 454, there is provided a 
predicted trajectory which will show how it is expected the 
plant will act and how the predictive model will predict. It 
is also noteworthy that the first value of the predicted 
trajectory is the actual input to the plant. The trajectory for 
the MV is computed every iteration using the previously 
described controller. The user therefore has the ability to 
view not only the actual response of the MV from a 
historical standpoint but, also the user can determine what 
the future prediction will be a number of steps into the 
future. 

[0280] A corresponding controlled variable (CV) curve is 
illustrated in FIG. 42. In this figure, there is also a historical 
section 450 and a predictive section 452. In the output, there 
is provided therein a desired response 462 which basically is 
a step response. In general, the system is set to vary from an 
output value of 80.0 to an input value of approximately 42.0 
which, due to dynamics, cannot be achieved in the real 
world. Various constraints arc also illustrated with a upper 
fuzzy constraint at a line 464 and a lower fuzzy constraint at 
a line 466. The system will show in the historical section the 
actual value on a bold line 468 which illustrates the actual 
response of the plant, this noted above as being either a 
simulator or by the plant itself. It should be remembered that 
the user can actually apply the input MV while the plant is 
running. At the current time line 454, a predicted value is 
shown along a curve 470. This response is what is predicted 
as a result of the input varying in accordance with the 



trajectory of FIG. 41. In addition to the upper and lower 
fuzzy constraints, there are also provided upper and lower 
hard constraints and upper and lower frustum values, that 
were described hereinabove. These are constraints on the 
trajectory which are determined during optimization. 

[0281] By allowing the user to view not only historical 
values but future predicted values during the operation of a 
plant or even during the simulation of a plant, the operator 
is now provided with information as to how the plant will 
operate during a desired change. If the operation falls 
outside of, for example, the upper and lower frustum of a 
desired operation, it is very easy for the operator to make a 
change in the desired value to customize this. This change is 
relatively easy to make and is made based upon future 
predicted behavior as compared to historical data. 

[0282] Referring now to FIG. 41, there is illustrated a 
block diagram of a control system for a plant that incorpo- 
rates the GUI interface. Basically, there is illustrated a plant 
480 which has a distributed control system (DCS) 482 
associated therewith for generating the MVs. The plant 
outputs the control variables (CV). A predictive controller 
484 is provided which is basically the controller as noted 
hereinabove with respect to FIG. 2 utilized in a control 
environment for predicting the future values of the manipu- 
lated variables MV(t+l) for input to the DCS 482. This will 
generate the predicted value for the next step. The predictive 
controller requires basically a number of inputs, the MVs, 
the output CVs and various other parameters for control 
thereof. A GUI interface 490 is provided which is operable 
to receive the MVs, the CVs, the predicted manipulated 
variables MV(t+l) for t+1, t+2, t+3, . t+n, from the predic- 
tive controller 484. It is also operable to receive a desired 
control variable CV D . The GUI interface 490 will also 
interface with a display 492 to display information for the 
user and allow the user to input information therein through 
a user input device 494, such as a mouse or pointing device. 
The predictive controller also provides to the GUI interface 
490 a predicted trajectory, which constitutes at one point 
thereof MV(t+l). The user input device 494, the display 492 
and the GUI interface 490 are generally portions of a 
software program that runs on a PC, as well as the predictive 
controller 484. 

[0283] In operation, the GUI interface 490 is operable to 
receive all of the information as noted above and to provide 
parameters on one or more lines 496 to the predictive 
controller 484 which will basically control the predictive 
controller 484, This can be in the form of varying the upper 
and lower constraints, the desired value and even possibly 
parameters of the identifying model. The model itself could, 
in this manner, be modified with the GUI interface 490. 

[0284] The GUI interface 490 allows the user to identify 
the model as noted hereinabove with respect to FIG. 39 and 
also allows the user to view the system in either simulation 
mode or in run time mode, as noted in FIG. 41. For the 
simulation mode, a predictive network 498 is provided 
which is operable to receive the values MV(t+l) and output 
a predicted control variable rather than the actual control 
variable. This is described hereinabove with reference to 
FIG. 2 wherein this network is utilized in primarily a 
predictive mode and is not utilized in a control mode. 
However, this predictive network must have the dynamics of 
the system model defined therein. 
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[0285] Referring dow to FIGS. 42-45, there is illustrated a 
screen view which comprises the layout screen for display- 
ing the variables, both input and output, during the simula- 
tion and/or run-time operation. In the view of FIG. 42, there 
is illustrated a set-up screen wherein four variables can be 
displayed: reflex, steam, top_comp, and bot_comp. These 
are displayed in a box 502. On the right side of the screen 
are displayed four boxes 504, 506, 508 and 510, for dis- 
playing the four variables. Each of the variables can be 
selected for which box they will be associated with. The 
boxes 504-510 basically comprise the final display screen 
during simulation. The number of rows and columns can be 
selected with boxes 512. 

[0286] Referring now to FIG. 42, there is illustrated the 
screen that will result when the system is "stepped" for the 
four variables as selected in FIG. 43. This will result in four 
screens displaying the two MVs, reflex and steam, and the 
two CVs, top_comp and bot_comp. In some situations, there 
may be a large number of variables that can be displayed on 
a single screen; there could be as many as thirty variables 
displayed in the simulation mode as a function of time. In the 
embodiment of FIG. 44, there are illustrated four simulation 
displays 516, 518, 520 and 522, associated with boxes 
504-510, respectively. 

[0287] Referring now to FIG. 43, there is illustrated 
another view of the screen of FIG. 42 with only one row 
selected with the boxes 512 in two columns. This will result 
in two boxes 524 and 528 that will be disposed on the final 
display. It can be seen that the content of these two boxes, 
after being defined by the boxes 512, is defined by moving 
to the variable box 502 and pointing to the appropriate one 
of the variable names, which will be associated with the 
display area 524 or 528. Once the variables have been 
associated with the particular display, then the user can 
move to the simulation screen illustrated in FIG. 45, which 
only has two boxes 530 and 532 associated with the boxes 
524 and 528 in FIG. 44. Therefore, the user can very easily 
go into the set-up screen of FIGS. 42 and 44 to define which 
variables will be displayed during the simulation process or 
during run-time execution. 

[0288] On-Line Optimizer 

[0289] Referring now to FIG. 46, there is illustrated a 
block diagram of a plant 600 utilizing an on-line optimizer. 
This plant is similar to the plant described hereinabove with 
reference to FIG. 6 in that the plant receives inputs u(t) 
which are comprised of two types of variables, manipulat- 
able variables (Mv) and disturbance variables (DV). The 
manipulatable variables are variables that can be controlled 
such as flow rates utilizing flow control valves, flame 
controls, etc. On the other hand the disturbance variables are 
inputs that are measurable but cannot be controlled such as 
feed rates received from another system. However, they all 
constitute part of the input vector u(t). The output of the 
plant constitutes the various states that are represented by 
the vector y(t). These outputs are input to an optimizer 602 
which is operable to receive desired values and associated 
constraints and generate optimized input desired values. 
Since this is a dynamic system, the output of the optimizer 
602 is then input to a controller 604 which generates the 
dynamic movements of the inputs which are input to the 
distributed control system 606 for generation of the inputs to 
the plant 600, it being understood that the DCS 606 will only 



generate the MVs that arc actually applied to the plant. The 
distinction of this system over other systems is that the 
optimizer 602 operates on-line. This aspect is distinctive 
from previous system in that the dynamics of the system 
must be accounted for during the operation of this optimizer. 
The reason for this is that when an input value moves from 
one value to another value, there are dynamics associated 
therewith. These dynamics must be considered or there will 
multiple errors, This is due to the fact that most predictive 
systems utilizing optimizers implement the optimization 
routine with steady state models. No decisions can be made 
until the plant settles out, such that such an optimization 
must be performed off-line. 

[0290] Referring now to FIG. 47, there is illustrated a 
block diagram of the optimizer 602. In general, there is 
provided a linear or nonlinear steady state model 610 which 
is operable to receive desired control variables (CV) which 
are basically the outputs of the plant y(t). The steady state 
model 610 is operable to receive the desired CV values as an 
input in the form of set points and generates an output 
optimized value of u(t+l). This represents the optimized or 
predicted future value that is to be input to the controller for 
controlling the system. However, the steady state model 610 
was generated with a training set of input vectors and output 
vectors that represented the plant at the time that this training 
data was taken. If, at a later time, the model became 
inaccurate due to changes in external uncontrollable aspects 
of the plant 600, then the model 610 would no longer be 
accurate. However, it is noted that the gain of the steady state 
model 610 would remain accurate due to the fact that an 
offset would be present. In order to account for this offset, 
a linear, non-linear dynamic model 612 is utilized which 
receives the inputs values and generates the output values 
CV which provide a prediction of the output value CV This 
is compared with an actual plant output value which is 
derived from an on-line analyzer or a virtual on-line ana- 
lyzer (VOA) 616. The VOA 616 is operable to receive the 
plant outputs CV and the plant states and the inputs u(t). The 
output of the VOA 616 provides the actual output of the 
plant which is input a difference circuit 618, the difference 
thereof being the offset or "bias." This bias represents an 
offset which is then filtered with a filter 620 for input to an 
offset device 622 to offset external CV set points, i.e., 
desired CV values, for input to the nonlinear steady state 
model 610. It should be noted that the VOA 616 could be the 
same steady state model 610. Also, the dynamic model 612 
can be replaced with an actual plant instrument reading. 

[0291] Referring now to FIG. 48, there is illustrated a 
diagrammatic view of the application of the bias illustrated 
in FIG. 47. There is illustrated a curve 630 representing the 
mapping of the input space to the output space through a 
steady state model. As such, for vector (MV), there will be 
provided a set of output variables (CV). Although the steady 
state model is illustrated as a straight line, it could have a 
more complex surface. Also, it should be understood that 
this is represented as a single dimension, but this could 
equally apply to a multi-variable system wherein there are 
multiple input variables for a given input vector and multiple 
output values for a given output vector. If the steady state 
model represents an accurate predictive model of the system, 
then input vector MV will correspond to a predicted output 
value. Alternatively, in a control environment, it is desirable 
to predict the MVs from a desired output value CV SET . 
However, the predicted output value or the predicted input 
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value, depending upon whether the input is predicted or the 
output is predicted, will be a function of the accuracy of the 
model. This can change due to various external unmeasur- 
able disturbances such as the external temperature, the 
buildup of slag in a boiler, etc. This will effectively change 
the way the plant operates and, therefore, the model will no 
longer be valid. However, the "gain" or sensitivity of the 
model should not change due to these external disturbances. 
As such, when an MV is varied, the output would be 
expected to vary from an initial starting point to a final 
resting point in a predictable manner. It is only the value of 
the predicted value at the starting point that is in question. 
In order to compensate for this, some type of bias must be 
determined and an offset provided. 

[0292] In the diagrammatic view of FIG. 48, there is 
provided a curve 632 representing the actual output of the 
plant as determined by the VOA 616 or from an actual plant 
instrument. A bias is measured at the MV point 634 such that 
an offset can then be determined. This offset will result in the 
steady state model being shifted upward across the model 
space by the bias value, as indicated by a dotted line 633. 
Therefore, for the CV SET value, the desired output value, the 
new MV value will be MV OPT . It can be seen that the gain 
of the model has not changed. 

[0293] Referring now to FIG. 49, there is illustrated a plot 
of both CV and Mv illustrating the dynamic operation. There 
is provided a first plot or curve 640 illustrating the operation 
of the CV output in response to a dynamic change which is 
the result of a step change in the MV input, represented by 
curve 642. The dynamic model 612 will effectively predict 
what will happen to the plant when it is not in steady state. 
This occurs in a region 646. In the upper region, a region 
648, this constitutes a steady state region. 

[0294] Referring now to FIG. 50, there is illustrated a plot 
illustrating the predicted model as being a steady state 
model. Again, there is provided the actual output curve 640 
which represents the output of the VOA or the value of the 
CV provided by an actual plant instrument and also a curve 
650 which represents the prediction provided by the steady 
state nonlinear model. It can be seen that, when the curve 
642 makes a step change, the output, as represented by the 
curve 640, will change gradually up to a steady state value. 
However, the steady state model will make an immediate 
calculation of what the steady state value should be, as 
represented by a transition 652. This will rise immediately to 
the steady state level such that, during the region 646, the 
prediction will be inaccurate. This is representative in the 
plot of bias for both the dynamic and the steady state 
configurations. In the bias for the steady state model, rep- 
resented by a solid curve 654, the steady state bias will 
become negative for a short period of time and then, during 
the region 646, go back to a bias equal to that of the dynamic 
model. The dynamic model bias is illustrated by a dotted line 
656. This example assumes a first order dynamic response, 
but higher order dynamic responses are equally applicable. 

[0295] Referring now to FIG. 51, there are illustrated 
curves representing the operation of the dynamic model 
wherein the dynamic model does not accurately predict the 
dynamics of the system. Again, the output of the system, 
represented by the VOA or an actual instrument reading of 
the CV, is represented by the curve 640. The dynamic model 
provides a predictive output incorporated in the dynamics of 



the system, which is represented by a curve 658. It can be 
seen that, during the dynamic portion of the curve repre- 
sented by region 646, that the dynamic model reaches a 
steady state value too quickly, i.e., it does not accurately 
model the dynamics of the plant during the transition in 
region 646. This is represented by a negative value in the 
bias, represented by curve 660 with a solid line. The filter 
620 is utilized to filter out the fast transitions represented by 
the dynamic model in the output bias value (not the output 
of the dynamic model itself). Also, it can be seen that the 
bias, represented by a dotted line 664, will be less negative. 

[0296] Referring now to FIG. 52, there is illustrated a 
block diagram of a prior art system utilizing a steady state 
model for optimization. In this system, there is provided a 
steady state model 670 similar to the steady state model 610. 
This is utilized to receive set points and predict optimized 
MVs, represented as the vector u OPT . However, in order to 
provide some type of bias for the operation thereof, the 
actual set points are input to an offset circuit 672 to be offset 
by a bias input. This bias input is generated by comparing 
the output of a steady state model 672, basically the same 
steady state model 670, with the output of a VOA 674, 
similar to VOA 616 in FIG. 47. This will provide a bias 
value, as it does in FIG. 47. However, it is noted that this is 
the bias between the steady state model, a steady state 
predicted value, and possibly the output of the plant which 
may be dynamic in nature. Therefore, it is not valid during 
dynamic changes of the system. It is only valid at a steady 
state condition. In order to utilize the bias output by a 
difference circuit 676 which compares the output of a steady 
state model 672 with that of the VOA 674, a steady state 
detector 678 is provided. This steady state detector 678 will 
look at the inputs and outputs of the network and determine 
when the outputs have "settled" to an acceptable level 
representative of a steady state condition. This can then be 
utilized to control a latch 680 which latches the output of the 
difference circuit 676, the bias value, which latched value is 
then input to the offset circuit 672. It can therefore be seen 
that this configuration can only be utilized in an off-line 
mode, i.e., when there are no dynamics being exhibited by 
the system. 

[0297] Steady State Optimization 

[0298] In general, steady state models are utilized for the 
steady state optimization. Steady state models represent the 
steady state mapping between inputs of the process (manipu- 
lated variables (MV) and disturbance variables (DV)) and 
outputs (controlled variables (CV)). Since the models rep- 
resent a steady state mapping, each input and output process 
is represented by a single input or output to the model (lime 
delays in the model are ignored for optimization purposes). 
In general, the gains of a steady state model must be accurate 
while the predictions are not required to be accurate. Pre- 
cision in the gain of the model is needed due to the fact that 
the steady state model is utilized in an optimization con- 
figuration. The steady state model need not yield an accurate 
prediction due to the fact that a precise VOA or actual 
installed instrument measurement of the CV can be used to 
properly bias the model. Therefore, the design of the steady 
state model is designed from the perspective of developing 
an accurate gain model during training thereof to capture the 
sensitivities of the plant. The model described hereinbelow 
for steady state optimization is a nonlinear model which is 
more desirable when modeling processes with multiple 
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operating regions. Moreover, when modeling a process with 
a single operating region, a linear model could be utilized. 
A single operating region process is defined as a process 
whose controlled variables operate at constant set-points, 
whose measured disturbance variables remain in the same 
region, and whose manipulated variables are only changed 
to reject unmeasured disturbances. Since the MVs are only 
moved to reject external disturbances, the process is defined 
to be external disturbance dominated. An example of a 
single operating region process is a distillation column 
operating at a single purity setpoint. By comparison, a 
multiple operating region process is a process whose 
manipulated variables are not only moved to reject unmea- 
sured disturbances, but are also changed to provide desired 
process performance. For example, the manipulated vari- 
ables may be changed to achieve different CV set points or 
they may be manipulated in response to large changes to 
measured disturbances. An example of this would be grade 
change control of a polypropylene reactor. Since the MVs or 
CVs of a multiple operating region process are often set to 
non-constant references to provide desired process perfor- 
mance, a multiple region process is reference dominated 
rather than disturbance dominated. The disclosed embodi- 
ment herein is directed toward multiple operating region 
processes and, therefore, a non- linear steady state model will 
be utilized. However, it should be understood that both the 
steady state model and the dynamic model could be linear in 
nature to account for single operating region processes. 

[0299] As described hereinabove, steady state models are 
generated or trained with the use of historical data or a first 
principals model, even a generic model of the system. The 
MVs should be utilized as inputs to the models but the states 
should not be used as inputs to the steady state models. 
Using states as inputs in the steady state models, i.e., the 
states being the outputs of the plant, produces models with 
accurate predictions but incorrect gains. For optimization 
purposes, as described hereinabove, the gain is the primary 
objective for the steady state optimization model. 

[0300] Referring now to FIG. 53, there is illustrated a 
block diagram of a model which is utilized to generate 
residual or computed disturbance variables (CDVs). A 
dynamic non-linear state model 690 which provides a model 
of the states of the plant, the measurable outputs of the plant, 
and the MVs and DVs. Therefore, this model is trained on 
the dynamics of the measurable outputs, the states, and the 
MVs and DVs as inputs. The predicted output, if accurate, 
should be identical to the actual output of the system. 
However, if there are some unmeasurable external distur- 
bances which affect the plant, then this prediction will be 
inaccurate due to the fact that the plant has changed over that 
which was originally modeled. Therefore, the actual state 
values, the measured outputs of the plant, are subtracted 
from the predicted states to provide a residual value in the 
form of the CDVs. Thereafter, the computed disturbances, 
the CDVs, are provided as an input to a non-linear steady 
state model 692, illustrated in FIG. 54, in addition to the 
MVs and DVs. This will provide a prediction of the CVs on 
the output thereof. Non-linear steady state model 692, 
described hereinabove, is created with historical data 
wherein the states are not used as inputs. However, the 
CDVs provide a correction for the external disturbances. 
This is generally referred to as a residual activation network 
which was disclosed in detail in U.S. Pat. No. 5,353,207 



issued Oct. 4, 1994 to J. Kceler, E. Hartman and B. Fergu- 
son, which patent is incorporated herein by reference. 

[0301] Referring now to FIG. 53a, there is illustrated a 
block diagram of the general residual activation neural 
network (RANN) architecture. As set forth in U.S. Pat. No. 
5,353,207, the state values are modeled in a state model 695 
which is comprised of an input layer, an output layer and a 
mapping layer wherein the input layer is operable to receive 
the independent values MV and DV. This corresponds to the 
state model 690. This provides from the output thereof 
predicted states. This state model generally provides a 
learned average of the state values. This is then input to a 
difference block 697 which compares the predicted value for 
the states with the actual values. It is important to note that 
the actual state values that are input to the block 697 only 
comprise those that are considered to be dependent upon the 
unmeasurable disturbance variables. For example, one of the 
unmeasurable disturbance variables is the quality of the coal 
input to the plant. From general engineering knowledge, it is 
known that steam temperature is one measurable state value 
having a value that exhibits a strong dependency upon the 
quality of the coal; that is, if the coal quality varies, then the 
steam temperature will vary. Although there may be empiri- 
cal relationships between the steam temperature and the coal 
quality, this is generally difficult to model. Therefore, if the 
coal quality were to vary from that predicted by the state 
model 695, this would result in a calculated difference and 
this calculated difference would constitute the unmeasurable 
disturbance which would be output by the difference block 
697 as the calculated disturbance variable (CDV). This is 
then input to a filter 699, which low pass filters out noise 
variations in the CDVs. The output of the filter 699, the 
filtered CDVs, is then input a main neural network, neural 
network 701, which receives both the MVs, DVs, and 
filtered CDVs. 

[0302] The use of the RANN in architecture in FIG. 53a 
allows the network model to incorporate some of the unmea- 
surable disturbances therein. This, therefore, allows the user 
to have a model that takes into account some of these 
unmeasurable concepts, which cannot be quantified directly 
from the data. As noted hereinabove, the concepts that are 
now incorporated into the model may include such things as 
variations in fuel quality and changes in boiler cleanliness 
(slag buildup) and the cleanliness of the Economizer and Air 
Preheater, as well as other key disturbance parameters for 
which there can be no direct on-line measurement. With the 
use of the slate model, an empirical model of the state 
variables is provided based on independent inputs and any 
variation between the original independent inputs and the 
model can be quantified with the difference block 697. 

[0303] On-Line Dynamic Optimization 

[0304] Referring now to FIG. 55, there is illustrated a 
detailed block diagram of the model described in FIG. 47. 
The non-linear steady state model 610, described herein- 
above, is trained utilizing manipulated variables (Mv) as 
inputs with the outputs being the CVs. In addition it also 
utilizes the DVs. The model therefore is a function of both 
the MVs, the DVs and also the CDVs, as follows: 

CV-f{MV f DV,CDV) 
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[0305] The noo-lincar steady state model 610 is utilized in 
an optimization mode wherein a cost function is defined to 
which the system is optimized such that the MVs can move 
to minimize a cost function while observing a set of con- 
straints. The cost function is defined as follows: 

J-ftMV, DV, CVset) 

[0306] noting that many other factors can be considered in 
the cost function, such as MV constraints, economic factors, 
etc. The optimized manipulatable variables (MV OPT ) are 
determined by iteratively moving the MVs based upon the 
derivative dJ/dMv. This is a conventional optimization tech- 
nique and is described in Mash, S. G. and Sofer, A., "Linear 
and Nonlinear Programming," McGraw Hill, 1996, which is 
incorporated herein by reference. This conventional steady 
state optimizer is represented by a block 700 which includes 
the non-linear steady state model 610 which receives both 
the CDVs, the DVs and a CV set point. However, the set 
point is offset by the offset block 672. This offset is deter- 
mined utilizing a non-linear dynamic prediction network 
comprised of the dynamic non-linear state model 690 for 
generating the CDVs, from FIG. 53, which CDVs are then 
input to a non-linear dynamic model 702. Therefore, the 
combination of the dynamic non-linear state model 690 for 
generating the CDVs and the non- linear dynamic model 702 
provide a dynamic prediction on output 704. This is input to 
the difference circuit 618 which provides the bias for input 
to the filter 620. Therefore, the output of the VOA 616 which 
receives both states as an input the MVs and DVs as inputs 
provides an output that represents the current output of the 
plant. This is compared to the predicted output and the 
difference thereof constitutes a bias. VOA 616 can be a real 
time analyzer that provides an accurate representation of the 
current output of the plant. The non-linear dynamic model 
702 is related to the non-linear steady state model 610 as 
described hereinabove, in that the gains are related. 

[0307] The use of the non- linear dynamic model 702 and 
the dynamic non-linear state model 690 provides a dynamic 
representation of the plant which can be compared to the 
output of the VOA 616. Therefore the bias will represent the 
steady state bias of the system and, therefore, can be utilized 
on line. 

[0308] Referring now to FIG. 56, there is illustrated a 
diagrammatic view of a furnace/boiler system which has 
associated therewith multiple levels of coal firing. The 
central portion of the furnace/boiler comprises a furnace 720 
which is associated with the boiler portion. The furnace 
portion 720 has associated therewith a plurality of delivery 
ports 722 spaced about the periphery of the boiler at different 
vertical elevations. Each of the delivery ports 722 has 
associated therewith a pulverizer 724 and a coal feeder 726. 
The coal feeder 726 is operable to feed coal into the 
pulverizer 724 at a predetermined rate. The pulverizer 724 
crushes the coal and mixes the pulverized coal with air. The 
air carries the coal to the delivery ports 722 to inject it into 
the furnace portion 720. The furnace/boiler will circulate the 
heated air through multiple boiler portions represented by a 
section 730 which provides various measured outputs (CV) 
associated with the boiler operation as will be described 
hereinbelow. In addition, the exhaust from the furnace which 
is recirculated, will have nitrous oxides (NOx) associated 
therewith. An NO x sensor 732 will be provided for that 
purpose. 



[0309] Referring now to FIG. 57 and FIG. 57a, there is 
illustrated a top cross-sectional view and a side cross- 
sectional view of the furnace portion 720 illustrating four of 
the delivery ports 722 spaced about the periphery of the 
furnace portion 720. The pulverized coal is directed into the 
furnace portion 720 to the interior thereof tangential to the 
furnace center point forming a swirling fireball interior to 
the furnace. This is what is referred to as a "tangentially 
fired" boiler. However, other boilers not utilizing tangential 
firing could be accommodated. Utilizing this technique, a 
conventional technique, a fireball can be created proximate 
to the delivery or inlet ports 722. Further, the feed rates for 
each of the elevations and the associated inlet ports 722 can 
be controlled in order to define the mass center of the fireball 
elevation. By varying the feed rates to the various eleva- 
tions, this fireball can be raised or lowered in the furnace. 
The placement of this fireball can have an effect on the 
efficiency, termed the "heat rate," the NO x level, and the 
output of unburned carbon (Loss on Ignition(LOI)). For this 
particular application, three important features to control or 
optimize are the heat rate, the NO x levels, and the LOI. 

[0310] In these boilers with multiple elevations of coil 
firing, the combination of the elevations in service and the 
amount of coal fed at each elevation are important param- 
eters when utilizing a prediction and optimization system. 
This is specifically so with respect to optimizing the NO x 
emissions and the performance parameters. Additionally, 
most of these boilers in the field have excess installed 
capacity for delivering fuel for each elevation and, therefore, 
opportunities exist to alter the method by which fuel is 
introduced using any given combination, as well as biasing 
the fuel up or down at any elevation. In general, for any 
given output level a relatively stable coal feed rate is 
required such that the increase or decrease of fuel flow to 
one elevation results in a corresponding, opposite direction 
change in coal flow to another elevation. A typical utility 
boiler will have between four to eight elevations of fuel 
firing and may have dual furnaces. This presents a problem 
in that representation of a plant in a neural network or some 
type of first principals model will require the model to 
represent the distribution of fuel throughout the boiler in an 
empirical model with between four and sixteen highly 
correlated, coal flow input variables. Neural networks, for 
example, being nonlinear in nature, will be more difficult to 
train with so many correlated variables. 

[0311] Much of the effect on the NO x emissions and 
performance parameters, due to these changes in fuel dis- 
tribution, relate to relative height in the boiler that the fuel 
is introduced. Concentrating the fuel in the bottom of the 
furnace by favoring the lower elevations of coal firing will 
yield different output results than that concentrating the fuel 
at the top of the furnace. The concept of Fuel Elevation has 
been developed in order to represent the relative elevation of 
the fuel in the furnace as a function of the feed rate and the 
elevation level. This provides a single parameter or manipu- 
latable variable for input to the network which is actually a 
function of multiple elevations and feed rate. The Fuel 
Elevation is defined as a dimensionless number that 
increases between "O" and "I" as fuel is introduced higher 
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in the furnace. Fuel Elevation is calculated according to the 
following equation: 

{KQjRj ) + {K 2 )(R 2 ) + (K 3 ){R A ) . ■ . -f (*„)( R n ) 

6 ~ R j + /? 2 + f? 3 . . . +• R n 

where: Fe - Calculated Fuel Elevation 

K\ ... K n = Elevation Constant for elevation 1 to n 

(described hereinbelow) 

R\...R n - Coal Feed Rate for elevation 1 to n 

[0312] Constants for each elevation are calculated to rep- 
resent the relative elevation at which the coal is introduced. 
For example, for a unit with four elevations of fuel firing, 
each elevation represents 25% of the combustion area. The 
Fuel Elevations constants for the lowest to the highest 
elevations are 0.125, 0375, 0.625 and 0.875. 

[0313] Referring now to FIG. 58, there is illustrated a 
block diagram of the boiler/furnace 720 connected in a 
feedback or control configuration utilizing an optimizer 740 
in the feedback loop. A controller 742 is provided, which 
controller 742 is operable to generate the various manipu- 
latable variables of the input 744. The manipulatable vari- 
ables, or MVs, are utilized to control the operation of the 
boiler. The boiler will provide multiple measurable outputs 
on an output 746 referred to as the CVs of the system or the 
variables to be controlled. In accordance with the disclosed 
embodiment, some of the outputs will be input to the 
optimizer 740 in addition to some of the inputs. 

[0314] There are some inputs that will be directly input to 
the optimizer 740, those represented by a vector input 748. 
However, there are a plurality of other inputs, represented by 
input vector 750, which are combined via a multiple MV- 
single M V algorithm 752 for input to the optimizer 740. This 
algorithm 752 is operable to reduce the number of inputs and 
utilize a representation of the relationship of the input values 
to some desired result associated with those inputs as a 
group. In the disclosed embodiment, this is the Fuel Eleva- 
tion. This, therefore, results in a single input on a line 754 
or a reduced set of inputs. 

[0315] The optimizer 740 is operable to receive a target 
CV on a vector input 756 and also various constraints on 
input 758. These constraints are utilized by the optimizer, as 
described hereinabove. This will provide a set of optimized 
MVs. Some of these MVs can be directly input to the 
controller, those that are of course correlated to the input 
vector 748. The input vector or MV corresponding to the 
vector input 754 will be passed through a single MV- 
multiple MV algorithm 760. This algorithm 760 is basically 
the inverse of the algorithm 752. In general, this will 
represent the above-noted Fuel Elevation calculation. How- 
ever, it should be recognized that the algorithm 752 could be 
represented by a neural network or some type of model, as 
could the algorithm 760, in additional to some type of 
empirical model. Therefore, the multiple inputs may be 
reduced to a lessor number of inputs or single input via some 
type of first principals algorithm or some type of predeter- 
mined relationship. In order to provide these inputs to the 
boiler, they must be processed through the inverse relation- 



ship. It is important to note, as described hereinabove, that 
the optimizer 740 will operate on-line, to account for dis- 
turbances. 

[0316] Referring now to FIG. 59, there is illustrated a 
block diagram of the training system for training the opti- 
mizer 740 as neural networks. A general model 770 is 
provided which is any type of trainable nonlinear model. 
These models are typically trained via some type of back- 
propagation technique such that a training database is 
required, this represented by training database 772. A train- 
ing system 774 operates the model 770 such that it is trained 
on the outputs and the inputs. Therefore, inputs are applied 
thereto with target outputs representing the plant. The 
weights are adjusted in the model through an iterative 
procedure until the error between the outputs and the inputs 
is minimized. This, again, is a conventional technique. 
Additionally, most models will be trained utilizing gain 
bounded training as described hereinabove. The gain 
bounded training is desirable since the gain of the model is 
important. The reason for this is that it is utilized in a 
feedback path for purposes of control. If the only purpose of 
this was prediction, then the gain would not be a sensitive 
aspect of the network modeling. 

[0317] In the disclosed embodiment, since there is defined 
a relationship between multiple inputs to a single or reduced 
set of inputs, it is necessary to train the model 770 with this 
relationship in place. Therefore, the algorithm 752 is 
required to reduce the plurality of inputs on vector 750 to a 
reduce set of inputs on the vector 754, which, in the 
disclosed embodiment, is a single value. This will constitute 
a single input or a reduced set of inputs that replace the 
multiple inputs on vector 750, which input represents some 
functional relationship between the inputs and some desired 
or observed behavior of the boiler. Thereafter, the remaining 
inputs are applied to the model 770 in the form of the vector 
748. Therefore, once the model 770 is trained, it is trained 
on the representation generated by the algorithm in the 
multiple MV-single MV algorithm 752. 

[0318] Referring now to FIGS. 60 and 61, there is illus- 
trated a more detailed block diagram of the embodiment of 
FIG. 55. The portion of the embodiment illustrated in FIG. 
59 is directed toward that necessary to generate the CV bias 
for biasing the set points. There are a plurality of measured 
inputs (MV) that are provided for the conventional boiler. 
These are the feeder speeds for each of the pulverizers, the 
Close Coupled Over-fired Air (CCOFA) value comprising 
various dampers that are adjusted on-line are: other inputs 
Separated Over-fired Air (SOFA) values which also are 
represented in terms of a percent open of select dampers, a 
tilt value, which defines the tilt of the inlet ports for injecting 
the fuel, the Wind Box to Furnace Pressure (WB/Fur) which 
is the differential pressure between the sandbox and the 
furnace, CCOFA tilt, which comprises the tilts on the Close 
Coupled Overtired Air Ports, plus additional inputs which 
are not transformed but are fed directly into the model. All 
are utilized to generate input variables for a network. With 
respect to outputs, the outputs will be in the form of the NOx 
values determined by a sensor: the dry gas loss, the main 
steam temperature, and the loss on ignition (LOI) for both 
reheat and superheat. 

[0319] Of the inputs, the feeder speeds are input to a Fuel 
Elevation algorithm block 790 which provides a single 
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output on an output 792 which is referred to as the Fuel 
Elevation concept, a single value. In addition, the multiple 
Auxiliary air Damper positions are input to an auxiliary air 
elevation algorithm block 794, which also provides a single 
value representing auxiliary air elevation, this not being 
described in detail herein, it being noted that this again is a 
single value representing a relationship between multiple 
inputs and a desired parameter of the model. The CCOFA 
values for each of the dampers provide a representation of a 
total CCOFA value and a fraction CCOFA value, and rep- 
resented by an algorithm block 796. This is also the case 
with the SOFA representation in a block 798 and also with 
the WB/Fur representation wherein a pseudo curve is uti- 
lized and a delta value is determined from that pseudo curve 
based upon the multiple inputs. This is represented by a 
block 800. Additionally, there is provided in a block 801 the 
tilt representation of the CCOFA tilts, which constitute the 
difference between the CCOFA tilt and the burner tilts. The 
block 801 receives on the input thereof the burner tilts and 
the CCOFA tilt and provides on the output thereof a Delta 
CCOFA tilt. 

[0320] The output of all of the blocks 790, 794, 796, 798, 
800 and 801 provide the MVs in addition to the Stack 0 2 
value on a line 802. These are all input to a state prediction 
model 804 similar to the model 690 in FIG. 55. This model 
also receives the disturbance variables (DV) for the system, 
these not being manipulatable inputs to the boiler. Hie state 
prediction then provides the predicted states to a difference 
block 806 for determining the difference between predicted 
states and the measured state variable output by the boiler 
720 on a line 808 which provides the measured states. The 
difference block 806 provides the raw Computed Distur- 
bance Variables (CDV) which are filtered in a filter 810. This 
essentially, for the boiler, will be due to the slag and the 
cleanliness of the boiler. Thereafter, the CD Vs, the MVs and 
the DVs are input to a model 812, which basically is the 
nonlinear model 702. The output of model 812 provides the 
estimated CVs which are then compared with the measured 
CVs in a difference block 814 to provide through a filter 816 
the CV bias values. 

[0321] Referring further to FIG. 61, the CV bias values, 
each CV representing a vector, output by the filter 816 is 
input to a steady state model optimization block 818, in 
addition to the MVs, DVs, and the CDVs. In addition, there 
are provided various targets in the form of set points, various 
optimization weightings and various constraints such as MV 
minimum and maximum constraints, in addition to various 
economic constraints for each CV. This is substantially 
identical to the steady state optimizer 700 illustrated in FIG. 
55. It should be understood that this is a mixed optimization 
technique employing both target based and economic based 
target functions, which defines a cost function which oper- 
ates on the targets, constraints, and the optimization weight- 
ings in the cost function relationship and then essentially 
utilizes the derivative of the cost function to determine the 
move of the MV. This is conventional and is described in 
Mash, S. G. and Sofer, A., "Linear and Nonlinear Program- 
ming," McGraw Hill, 1996, which was incorporated herein 
by reference. 

[0322] The Auxiliary Air Elevation block 794 provides a 
relationship between the auxiliary air dampers and the 
concept of both providing a simple geometric center and 
flowated geometric center. As an example, consider that 



there arc provided three dampers as inputs, Damper 3 (D3), 
Damper 5 (D5), and Damper 7 (D7). Of course, additional 
dampers could be provided. With these three dampers, the 
following relationship will exist: 



{j^j(D3) + {D5) + {\){D7)} 



[0323] wherein each value for the damper represents a 
percent open value therefor. This represents the center of 
mass for the auxiliary air based on the damper positions and 
an arbitrary furnace elevation scale wherein in Damper 7 is 
set at 1, Damper 5 is set at 2 /i of the elevation, and Damper 
3 is set at Vi of the elevation. 

[0324] The flowrated geometric center concept is also 
associated with the auxiliary air center and is set forth as 
follows: 



Aux_air_center = - 

D 

where: 

a = ( l - j( D3){D3j\rea)[^ WB / FUR ) + 

( 2 - ^{D5){D5j\ rea)[ V WB / FUR ) + 

( 1 ){D7){D7_Area)( V WB/ FUR ) 
b = (D3){D3_A rea){y/ WB / FUR ) + 
(D5)(D5J{ rea){ V WB / FUR ) + 
(D7)(D7 'j\rea)(yf WB/FUR ) 

D3=% open of Damper 3 
D5=% open of Damper 5 
D7 = % open of Damper 7 
WB/FUR = Wind box to Furnace pressure - 

differential pressure between the windbox and 

the interior of the furnace 

[0325] The optimization model 818 will provide MV set 
points. These MV set points could be, for such MVs as the 
Stack 0 2 , input directly to the boiler 720 for control thereof 
as a new input value. However, when the MVs that represent 
the single values such as Fuel Elevation which relates back 
to multiple inputs must be processed through the inverse of 
that relationship to generate the multiple inputs. For 
example, Fuel Elevation values are provided as M V on a line 
820 for input to a Fuel Elevation neural network 822 which 
models the relationship of Fuel Elevation to feeder speed. 
However, a neural network is not necessarily required in that 
the basic relationship described hereinabove with respect to 
Fuel Elevation will, in some cases, suffice and the algorithm 
required is only the inverse of that relationship. In other 
cases, the feasible feeder combinations must be considered 
and a neural network with binary optimization employed. 
Binary optimization is an iterative technique whereby less 
than all of the feeders arc selected and that particular 
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combination analyzed. If, for example, there were n feeders, 
this would result in two 2 n combinations. Not all of these 2" 
combinations would be utilized. Only select ones would be 
utilized and these analyzed and optimized. After optimiza- 
tion of each of the selected combinations, they are analyzed 
and then one combination selected. In any of the discussed 
cases, this will provide on the output feeder speeds on lines 
824, which are multiple inputs. In addition, the auxiliary 
elevation is processed through a representation of a neural 
network which relates the multiple auxiliary air open values 
to the M V input in a block 826. The CCOFA representation 
provides the same inverse relationship in a block 828 to 
provide the MV set points associated therewith, the CCOFA, 
and the fraction CCOFA and provide the CCOFA percent 
open values. Of the SOFA MVs, the total SOFA and the 
fraction SOFA are processed through an inverse SOFA 
representation to provide the SOFA percent open inputs to 
the boiler 720. Lastly, the delta value from the WB/Fur curve 
is provided as MV set point through an inverse relationship 
in a block 832 in order to determine the WB/Fur input value 
to the boiler 720. The Delta CCOFA tilt value is processed 
through an inverse relationship of the tilt representation in a 
block 831 to provide on the output thereof the CCOPA tilt 
value. All of these operations, the optimization and the 
conversion operations, are done in real time, such that they 
must take into account the real-time variations of this 
system. Further, as described hereinabove, by reducing the 
amount of inputs, the actual steady state models and 
dynamic models will provide a better representation, and the 
sensitivities have been noted as being augmented for these 
representations. With this technique, the center of mass of 
the ball in the furnace 720 can be positioned with the use of 
one representative input modeled in the neural network or 
similar type model to allow efficiency and NOx to be 
optimized. It is noted that each of the inputs that represents 
multiple inputs to any of the algorithm blocks noted here- 
inabove and which are represented by a single variable each 
have a predetermined relationship to each other, i.e., the feed 
rate at each elevation has some relationship to the other 
elevations only with respect to a parameter defined as the 
center of mass of the fireball. Otherwise, each of these feeder 
rates is independent of each other. By defining a single 
parameter that is of interest and the relationship between the 
inputs to define this relationship, then that parameter itself 
can become a more important aspect of the model. 

[0326] Referring now to FIG. 62, there is illustrated a 
diagrammatic view of a process depicting the concept of 
taking multiple input values and converting them through a 
mapping relationship to a single value, processing that value 
to provide a modified value and then expanding that single 
modified value or concept back into the individual inputs as 
modified values therefore in accordance with the process. 
The single concept, as described hereinabove, is referred to 
as Fe. The inputs are referred to as values x 1( x 2 . . . x n 
referring to inputs that compromise more than one input. 
This procedure can be applied to as few as two inputs, three 
inputs, or any number of inputs, the feature being that 
multiple inputs are reduced per some relationship to less 
than the number of multiple inputs. These multiple inputs 
are processed through a mapping layer 900 to provide a 
single output Fe. It is important to note that the mapping 
process 900 is operable to reduce the multiple inputs to a 
number less than the number of multiple inputs. Therefore, 
if there were N inputs, and O outputs, the constraint is that 



O is less than N. There could be ten inputs with a single 
output or ten inputs with nine outputs. The important fact is 
that there is some relationship in the mapping layer 900 
which relates the output or outputs back to the inputs. 

[0327] Once the inputs have been mapped through the 
mapping layer 900, then the output is subjected to some type 
of process through a process block 902. In the disclosed 
embodiment, this is an optimization step. However, any type 
of process can be performed on the output of the mapping 
layer 900 to modify the value thereof, it being noted that 
operating on less outputs than the total number of inputs 
provides a different level of latitude to the model designer. 
This will provide on the output thereof a modified output Fe'. 
This is then expanded to an inverse mapping layer 904 
which then provides the multiple modified input values x/, 
x 2 ' • • • x n \ Each of these values x/, x^ . , . is then 
processed through individual process blocks 906. Individual 
process blocks 906 are process blocks that are associated 
only with that particular modified input. This is to be 
compared to a situation wherein the process defined in the 
process block 902 and the process blocks 906 would be 
combined and applied individually to each one of the inputs 
Xj, x 2 . . . x n . With the use of the system disclosed herein 
wherein a single concept can be defined by some relation- 
ship dependant upon multiple inputs, this concept can be 
processed by itself without regard to the individual inputs 
(other than the inherent depend ancies of the process on the 
inputs). Thereafter, once the desired behavior of the concept 
is achieved through the process block 902, i.e., through an 
optimization step, then the individual inputs can be pro- 
cessed separately, if necessary, 

[0328] Referring now to FIG. 63, there is illustrated an 
application of the diagrammatic view of FIG. 62. In this 
figure, there is illustrated a block 920 which represents the 
mapping layer 900 for receiving the inputs x 1( x 2 . . . x n . The 
mapping block 920 maps each of the inputs x u . . . x„ 
through the relationship described hereinabove for fuel 
elevation as follows: 

Fe^F{x v x^ . . . 

[0329] This provides on the output the single output F e . 
Again, it should be noted that there could be more than one 
output, depending upon the relationship map in the block 
920. 

[0330] The single value Fe representing the relationship in 
block 920 is then input to an optimizer block 922 which is 
operable to optimize the value according to various opti- 
mizer constraints. For example, if the fuel elevation were 
found to be slightly off from an optimization standpoint and 
it was desirable to calculate a new fuel elevation, then a new 
value would be output. Since the fuel elevation takes various 
feeder values and relates those to some type of spatial 
relationship of a fireball, this spatial relationship must then 
be converted back to feeder values. This is accomplished by 
processing this optimized value through the inverse rela- 
tionship represented in block 924 which provides a relation- 
ship: 

F~ l (x v x2, . . . xj 

[0331] This provides the outputs x^', x 2 . . . x^'. These 
represent the feeder values that would provide the new fuel 
elevation position. Of course, this is a predicted value and 
may not be correct. Of course, during the optimization 
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operation, this might require a number of passes through the 
optimization block 922. However, there are certain situa- 
tions wherein constraints and limits must be placed upon the 
outputs to ensure that they are within the safe operating 
limits of the control system. In this manner, post -optimizing 
constraints are placed on the outputs x x \ x 2 . . . x n ' in a block 
926. This represents the process blocks 906 in FIG. 62. This 
provides on the output x/, x 2 " , . . x n ". These constraints are 
such things as hard limits, soft limits, fuzzy limits, combi- 
natorial constraints, etc., which were described in U.S. Pat. 
No. 5,781,432, issued Jul. 14, 1998 to the present Assignee, 
which patent is incorporated herein by reference. 

[0332] In order to adequately deal with a single conceptual 
value, Fe, it is necessary to ensure that this value, during 
optimization, does not result in an "unobtainable" value, due 
to the inability of the inputs to operate over a wide enough 
range to provide the optimized value of Fe. In order to 
ensure that the optimized value of Fe is within acceptable 
limits, constraints must be defined therefore. Since Fe 
defines a relationship among multiple values to a single 
value, it is necessary to define the constraints for Fe based 
upon the constraints associated with the various input values 
x lf x 2 - . . x n . The relationship F(x 1 , x 2 . . . x n ) is a linear 
relationship in the disclosed embodiment and constraints on 
any of the inputs x lf x 2 . . . x n can be accumulated and 
directly related to constraints on the output Fe. Essentially, 
for the relationships disclosed herein, each of the inputs is 
evaluated at their maximum constraint and their minimum 
constraint which thereby defines a maximum or minimum 
constraint on Fe. This is represented by a block 928 which 
then provides the optimizer constraints in a block 930. 
Therefore, it is now known that, if the optimization opera- 
tion in block 922 optimizes Fe to provide the optimized 
value Fe' value, that the Fe' value will be within the defined 
constraints acceptable for the inputs and, therefore, the 
values of Xj, x 2 . . . x n can be realized within those 
constraints. However, it is important that when the inverse 
relationship F" 1 ^, x 2 . . . xj is evaluated, it takes into 
account that these constraints are in place. This will result in 
the values Xj', x 2 , . . . x 2 \ . . . x n ' being within the constraint 
values as were defined by block 928. The post optimizing 
constraints in block 926 provide additional constraints that 
were not embodied within the constraints defined in block 
920 or the embedded constraints in block 924. 

[0333] In summary, there has been provided a method for 
processing multiple inputs that can be correlated to desired 
performance parameters of a system. These performance 
parameters are defined in a relationship as a function of the 
multiple inputs. The inputs are first mapped through this 
relationship to provide intermediate values that number less 
than the number of inputs. These intermediate values are 
then modified in accordance with a predetermined modifi- 
cation algorithm to provide modified parameters. These 
modified parameters are then converted back to the input 
values through the inverse relationship such that these 
constitute modified inputs modified by the modification 
process. 

[0334] Although the preferred embodiment has been 
described in detail, it should be understood that various 
changes, substitutions and alterations can be made therein 
without departing from the spirit and scope of the invention 
as defined by the appended claims. 



What is claimed is: 

1. A method for modifying the values of select ones of a 
plurality of inputs to a plant, comprised of steps of: 

receiving the plurality of inputs; 

mapping in a first step of mapping N select ones of the 
received inputs through a first relationship defining P 
intermediate inputs numbering less than N, wherein the 
first relationship defines the relationship between the N 
select received inputs and the P intermediate inputs as 
a set of P intermediate values; 

processing the P intermediate values through a modifying 
process which modifies the P intermediate values in 
accordance with a predetermined modification algo- 
rithm to provide P modified intermediate inputs with 
corresponding P modified intermediate values; and 

mapping in a second step of mapping the P modified 
intermediate inputs and associated P modified interme- 
diate values through the inverse of the first relationship 
to provide N modified select input values at the output 
of the second mapping step corresponding to the N 
select received input values. 

2. The method of claim 1, and further comprising the step 
of post processing select ones of the N modified select input 
values on the N modified select inputs with a predetermined 
process algorithm. 

3. The method of claim 2, wherein the predetermined 
process algorithm is unique to the associated one of the N 
modified select inputs. 

4. The method of claim 2, wherein the predetermined 
process algorithms provide constraining limits on the range 
of associated ones of the N modified select inputs. 

5. The method of claim 1 , wherein the step of processing 
the P intermediate inputs through a modifying process 
comprises the step of processing the P intermediate input 
values through an optimizer to provide optimized interme- 
diate input values as the P modified intermediate values. 

6. The method of claim 5, wherein the optimizer is a 
steady state optimizer. 

7. The method of claim 6, wherein the steady state 
optimizer comprises a neural network which is non-linear. 

8. The method of claim 1, wherein the P intermediate 
inputs constitute a spatial relationship for a parameter of the 
plant wherein the first relationship defines the spatial rela- 
tionship of the P intermediate inputs to the N select received 
inputs. 

9. The method of claim 1 and further comprising the step 
of pre-processing select ones of the N select ones of the 
received inputs through a predetermined pre-process algo- 
rithm. 

10. The method of claim 9, wherein the step of pre- 
processing comprises: 

defining constraints for select ones of the N select ones of 
the received inputs; and 

limiting the predetermined modification algorithms such 
that processing of the P intermediate values will not 
result in the corresponding ones of the N modified 
inputs from exceeding the defined constraints. 
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11 . A method for optimization of the operation of a boiler, 
comprising the steps of: 

measuring the inputs and the outputs of the boiler; 

mapping a defined plurality of the measured inputs 
through a predetermined relationship that defines a 
desired operating parameter of the plant based upon 
said defined plurality of the measured inputs to inter- 
mediate inputs numbering less than the defined plural- 
ity of the measured inputs; 

processing the intermediate inputs and the inputs not in 
said defined plurality of the measured inputs through an 
optimizer to provide optimized intermediate input val- 
ues for the intermediate inputs and optimized inputs not 
in the defined plurality of the measured inputs opti- 
mized in accordance with a predetermined optimization 
algorithm; 

mapping the optimized intermediate input values through 
an inverse of the predetermined relationship to provide 
an optimized defined plurality of inputs corresponding 
to the defined plurality of the measured inputs; and 

applying the optimized defined plurality of inputs and the 
optimized inputs not in the defined plurality of the 
measured inputs to the boiler. 

12. The method of claim 11, wherein the boiler has a 
furnace with a plurality of inlet ports for providing fuel to a 
desired location in the furnace and arranged at different 
elevations, such that fuel can selectively be delivered at 
various of the elevations at different feed rates, and wherein 
the predetermined relationship comprises: 



) + (*l)(*2) + (*3)(*4) - . . + (K n )(R n ) 
/?| + R 2 +/?3 ... +/? n 



where: 

Fe=CaIculated Fuel Elevation 

Kj . . . K^Elevation Constant for elevation 1 to n 
(described hereinbelow) 

R a . . . R n =Coal Feed Rate for elevation 1 to n such that 
fuel can be concentrated at a desired location in the 
furnace with a desired center of mass for the fuel 
concentration. 

13. The method of claim 11, wherein the boiler has a 
furnace with a plurality of dampers for providing auxiliary 
air to the furnace with a geometric center defined therefor 
and arranged at different elevations, such that air can selec- 
tively be delivered at various of the elevations at different 
rates, and wherein the predetermined relationship com- 
prises: 

f(i](/>3) + U>5) + (l)(£>7)) 
Aux_air_center = D3 + D5+D7 

wherein D3, D5 and D7 represent damper opening values for 
three dampers, each value for the damper representing a 
percent open value therefor, this relationship representing 
the center of mass for the auxiliary air based on the damper 



positions and an arbitrary furnace elevation scale wherein in 
damper D7 is set at the maximum elevation, damper D5 is 
set at ft of the elevation, and damper D3 is set at x h of the 
elevation. 

14. The method of claim 11, wherein the boiler has a 
furnace with a plurality of dampers for providing auxiliary 
air to the furnace from a windbox with a geometric center 
defined therefor and arranged at different elevations, such 
that air can selectively be delivered at various of the eleva- 
tions at different rates, and wherein the predetermined 
relationship comprises: 



Aux_air_center = - 



where D3, D5 and D7 represent damper opening values for 
three dampers, each value for the damper representing a 
percent open value therefor, this relationship representing 
the flowated geometric center of for the auxiliary air and: 




\{D3){D3jirea){-J WB/FUR ) + 



^{D5)(D5j\rea){^ WB/FUR ) + 

( 1 ){D7){D7j\ rea)( V WB / FUR ) 
b = (D3){D3_Area){^ WB/FUR ) + 
(D5){D5_Area)(yf WB/ FUR ) + 
(D7){D7_Area)(^WB/ FUR ) 



D3% open of Damper 3 
D5«% open of Damper 5 
D7=% open of Damper 7 

WB/FUR=Windbox to Furnace pressure-differential pres- 
sure between the windbox and the interior of the 
furnace. 

15. A method for optimizing the operation of a plant 
having a plurality of manipulatable input variables (MVs) 
and at least one output, comprising the steps of: 

providing a predictive network to receive select ones of 
the MVs and predict new and updated values therefor 
as updated MV inputs to the plant in accordance with 
a predetermined prediction algorithm; and 

inputting the updated MV inputs to the plant; 

the step of providing a predictive network, comprising the 
steps of: 

mapping in a first step of mapping N select ones of the 
MVs through a first relationship defining P interme- 
diate MVs numbering less than N, wherein the first 
relationship defines the relationship between the N 
select MVs and the P intermediate MVs as a set of 
P intermediate values, 

processing the P intermediate values and the non 
mapped ones of the MVs through a modifying pro- 
cess which modifies the P intermediate values and 



06/03/2004, EAST Version: 1.4.1 



US 2003/0018399 Al Jan. 23, 2003 

31 



the non mapped ones of the MVs in accordance with 
a predetermined modification algorithm to provide P 
modified intermediate MVs with corresponding P 
modified intermediate values and non mapped modi- 
fied MVs for the non mapped ones of the MVs, and 

mapping in a second step of mapping the P modified 
intermediate MVs and associated P modified inter- 
mediate values through the inverse of the first rela- 
tionship to provide N modified select MVs at the 
output of the second mapping step corresponding to 
the N select MVs, 

wherein the P modified intermediate MVs and the non 
mapped modified MVs comprise the updated MV 
inputs. 

16. The method of claim 15, and further comprising the 
step of post processing select ones of the N modified MVs 
on the N modified MVs with a predetermined process 
algorithm. 

17. The method of claim 16, wherein the predetermined 
process algorithm is unique to the associated one of the N 
modified MVs. 

18. The method of claim 16, wherein the predetermined 
process algorithms provide constraining limits on the range 
of associated ones of the N modified select inputs. 

19. The method of claim 15, wherein the step of process- 
ing the P intermediate MVs and the non mapped MVs 
through a modifying process comprises the step of process- 
ing the P intermediate input values and the non mapped MVs 
through an optimizer to provide optimized intermediate 
MVs as the P modified intermediate values and optimized 
non mapped Mvs for the non mapped MVs. 

20. The method of claim 19, wherein the optimizer is a 
steady state optimizer. 

21. The method of claim 20, wherein the steady state 
optimizer comprises a neural network which is non-linear. 

22. The method of claim 15, wherein the P intermediate 
MVs constitute a spatial relationship for a parameter of the 
plant wherein the first relationship defines the spatial rela- 
tionship of the P intermediate MVs to the N select MVs. 

23. The method of claim 15 and further comprising the 
step of pre-processing select ones of the N select ones of the 
MVs through a predetermined pre-process algorithm. 

24. The method of claim 23, wherein the step of pre- 
processing comprises: 

defining constraints for select ones of the N select ones of 
the MVs; and 

limiting the predetermined modification algorithms such 
that processing of the P intermediate values will not 
result in the corresponding ones of the N modified MVs 
from exceeding the defined constraints. 

25. The method of claim 15, wherein there are unmea- 
surable disturbances associated with the plant and the modi- 
fying process is based on a model of the plant, and further 
comprising the step of determining an estimation of the 
unmeasurable disturbances and altering the modifying pro- 
cess based upon such estimation to account for the unmea- 
surable disturbances. 

26. The method of claim 25, wherein the step of process- 
ing the P intermediate MVs and the non mapped MVs 
through a modifying process comprises the step of process- 
ing the P intermediate input values and the non mapped MVs 
through a steady state optimizer to provide optimized inter- 



mediate MVs as the P modified intermediate values and 
optimized non mapped Mvs for the non mapped Mvs and the 
step of determining and altering operable to determine a bias 
value for the steady state optimizer and apply it to the steady 
state optimizer, 

27. A mapping network for modifying the values of select 
ones of a plurality of inputs to a plant, comprising: 

an input layer for receiving the plurality of inputs; 

an input mapper for mapping N select ones of said 
received inputs through a first relationship defining P 
intermediate inputs numbering less than N, wherein 
said first relationship defines the relationship between 
said N select received inputs and said P intermediate 
inputs as a set of P intermediate values; 

a processor for processing said P intermediate values 
through a modifying process which modifies said P 
intermediate values in accordance with a predeter- 
mined modification algorithm to provide P modified 
intermediate inputs with corresponding P modified 
intermediate values; and 

an output mapper for mapping said P modified interme- 
diate inputs and associated P modified intermediate 
values through the inverse of said first relationship to 
provide N modified select input values at the output of 
said output mapper corresponding to said N select 
received input values. 

28. The mapping network of claim 27, and further com- 
prising a post processor for processing select ones of said N 
modified select input values on said N modified select inputs 
with a predetermined process algorithm. 

29. The mapping network of claim 28, wherein the 
predetermined process algorithm is unique to the associated 
one of said N modified select inputs. 

30. The mapping network of claim 28, wherein the 
predetermined process algorithm provide constraining limits 
on the range of associated ones of said N modified select 
inputs. 

31. The mapping network of claim 27, wherein said 
processor is operable to process said P intermediate input 
values through an optimizer to provide optimized interme- 
diate input values as said P modified intermediate values. 

32. The mapping network of claim 30, wherein said 
optimizer is a steady state optimizer. 

33. The mapping network of claim 32, wherein said 
steady state optimizer comprises a neural network which is 
non-linear. 

34. The mapping network of claim 27, wherein said P 
intermediate inputs constitute a spatial relationship for a 
parameter of the plant wherein said first relationship defines 
the spatial relationship of said P intermediate inputs to said 
N select received inputs. 

35. The mapping network of claim 27 and further com- 
prising a pre-processor for processing select ones of said N 
select ones of said received inputs through a predetermined 
pre-process algorithm. 

36. The mapping network of claim 35, wherein said 
pre-processor comprises: 

a constraining device for defining constraints for select 
ones of said N select ones of said received inputs; and 

said constraining device for limiting said predetermined 
modification algorithms such that processing of said P 
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intermediate values by said processor will not result in 
the corresponding ones of said N modified inputs from 
exceeding said defined constraints. 

37. An optimization network for optimizing the operation 
of a boiler, comprising: 

an input network for measuring the inputs and the outputs 
of the boiler; 

an input mapper for mapping a defined plurality of said 
measured inputs through a predetermined relationship 
that defines a desired operating parameter of the plant 
based upon said defined plurality of said measured 
inputs to intermediate inputs numbering less than said 
defined plurality of said measured inputs; 

a processor for processing said intermediate inputs and 
said inputs not in said defined plurality of said mea- 
sured inputs through an optimizer to provide optimized 
intermediate input values for said intermediate inputs 
and optimized inputs not in said defined plurality of 
said measured inputs optimized in accordance with a 
predetermined optimization algorithm; 

an output mapper for mapping said optimized intermedi- 
ate input values through an inverse of said predeter- 
mined relationship to provide an optimized defined 
plurality of inputs corresponding to the defined plural- 
ity of said measured inputs; and 

a controller for applying said optimized defined plurality 
of inputs and said optimized inputs not in said defined 
plurality of said measured inputs to the boiler. 

38. The optimization network of claim 37, wherein the 
boiler has a furnace with a plurality of inlet ports for 
providing fuel to a desired location in said furnace and 
arranged at different elevations, such that fuel can selec- 
tively be delivered at various of the elevations at different 
feed rates, and wherein said predetermined relationship 
comprises: 

Fc ,. + iK 2 m) + (*3)(*4) - + (XJtfl.) 

6 /?| + R 2 ... + R„ 



where: 

Fc=CaIculatcd Fuel Elevation 

K n . . . R n =Elcvation Constant for elevation 1 to n 
(described hereinbelow) 

Rj . . . R n =Coal Feed Rate for elevation 1 to n such that 
fuel can be concentrated at a desired location in said 
furnace with a desired center of mass for said fuel 
concentration. 

39. The optimization network of claim 37, wherein the - 
boiler has a furnace with a plurality of dampers for providing 
auxiliary air to said furnace with a geometric center defined 
therefor and arranged at different elevations, such that air 
can selectively be delivered at various of t^e elevations at 
different rates, and wherein said predetermined relationship 
comprises: 



((i](0?) + (D5) + (l)(D7)] 

Aux_air_ccnter = — — — — — 

D3 + D5 + D7 

wherein D3, D5 and D7 represent damper opening values for 
three dampers, each value for said damper representing a 
percent open value therefor, this relationship representing 
the center of mass for the auxiliary air based on said damper 
positions and an arbitrary furnace elevation scale wherein in 
damper D7 is set at said maximum elevation, damper D5 is 
set at 2 A of the elevation, and damper D3 is set at Vz of the 
elevation. 

40. The method of claim 36, wherein the boiler has a 
furnace with a plurality of dampers for providing auxiliary 
air to said furnace from a windbox with a geometric center 
defined therefor and arranged at different elevations, such 
that air can selectively be delivered at various of the eleva- 
tions at different rates, and wherein said predetermined 
relationship comprises: 



Aux_ait_center= - 



where D3, D5 and D7 represent damper opening values for 
three dampers, each value for said damper representing a 
percent open value therefor, this relationship representing 
the flowated geometric center of for the auxiliary air and: 



a = ](D3)(D3 j\rea)[JWB/FUR) + 

( ? }{D5)(D5_Area)( V WB / FUR ) +• ( 1 ){D7){D7 _Area)(-J WB/ FUR ) 

b = {D3){D3j\rea)[y} WB/ FUR ) + {D5){D5j\rea){^ WB/ FUR ) + 

(D7)(D7 _Area)ftWB/FUR) 



D3~% open of Damper 3 
D5=% open of Damper 5 
D7«% open of Damper 7 

WB/FUR -Windbox to Furnace pressure — differential 
pressure between said windbox and the interior of said 
furnace. 

41. A optimization network for optimizing the operation 
of a plant having a plurality of manipulatable input variables 
(MVs) and at least one output, comprising; 

a predictive network operable to receive select ones of the 
MVs on an input thereto and predict and output on an 
output new and updated values therefor as updated M V 
inputs to the plant in accordance with a predetermined 
prediction algorithm; and 

a controller for inputting the updated MV inputs to the 
plant; 

said predictive network having: 

an input mapper for mapping N select ones of the MVs 
through a first relationship defining P intermediate 
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MVs numbering less than N, wherein said first 
relationship defines the relationship between said N 
select MVs and said P intermediate MVs as a set of 
P intermediate values, 

a processor for processing said P intermediate values 
and the non mapped ones of said MVs through a 
modifying process which modifies said P intermedi- 
ate values and said non mapped ones of said MVs in 
accordance with a predetermined modification algo- 
rithm to provide P modified intermediate MVs with 
corresponding P modified intermediate values and 
non mapped modified MVs for said non mapped 
ones of said MVs, and 

an output mapper for mapping said P modified inter- 
mediate MVs and associated P modified intermediate 
values through the inverse of said first relationship to 
provide N modified select MVs at the output of the 
second mapper corresponding to said N select MVs, 

wherein said P modified intermediate MVs and said 
non mapped modified MVs comprise said updated 
MV inputs. 

42. The optimization network of claim 40, and further 
comprising a post processor for processing select ones of 
said N modified MVs on said N modified MVs with a 
predetermined process algorithm. 

43. The optimization network of claim 41, wherein said 
predetermined process algorithm is unique to said associated 
one of said N modified MVs. 

44. The optimization network of claim 41, wherein said 
predetermined process algorithm provides constraining lim- 
its on the range of associated ones of said N modified select 
inputs. 

45. The method of claim 40, wherein said processor is 
operable to process said P intermediate input values and said 
non mapped MVs through an optimizer to provide optimized 
intermediate MVs as said P modified intermediate values 
and optimized non mapped Mvs for said non mapped MVs. 

46. The optimization network of claim 44, wherein said 
optimizer is a steady state optimizer. 



47. The optimization network of claim 19, wherein said 
steady state optimizer comprises a neural network which is 
non-linear. 

48. The meth optimization network of claim 40, wherein 
said P intermediate MVs constitute a spatial relationship for 
a parameter of the plant wherein said first relationship 
defines the spatial relationship of said P intermediate MVs 
to said N select MVs. 

49. The optimization network of claim 40 and further 
comprising a pre-processor for processing select ones of said 
N select ones of said MVs through a predetermined pre- 
process algorithm. 

50. The optimization network of claim 48, wherein the 
step of pre-processing comprises: 

defining constraints for select ones of said N select ones 
of said MVs; and 

limiting said predetermined modification algorithms such 
that processing of said P intermediate values will not 
result in the corresponding ones of said N modified 
MVs from exceeding said defined constraints. 

51. The optimization network of claim 40, wherein there 
are unmeasurable disturbances associated with the plant and 
said modifying process is based on a model of the plant, and 
further comprising an estimation network for determining an 
estimation of the unmeasurable disturbances and altering 
said modifying process based upon such estimation to 
account for said unmeasurable disturbances. 

52. The optimization network of claim 50, wherein said 
processor is operable to process said P intermediate input 
values and the non mapped MVs through a steady state 
optimizer to provide optimized intermediate MVs as said P 
modified intermediate values and optimized non mapped 
Mvs for said non mapped Mvs and said estimation network 
is operable to determine a bias value for said steady state 
optimizer and apply it to said steady state optimizer. 

* * * * * 
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